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ABSTRACT 

Using hydrodynamic simulations of disk-galaxy major mergers we investigate the star- 
formation history and remnant properties when various parameterizations of a simple 
stellar feedback model are implemented. The simulations include radiative cooling, a 
density-dependent star-formation recipe and a model for feedback from massive stars. 
The feedback model stores supernova feedback energy within individual gas parti¬ 
cles and dissipates this energy on a time-scale specified by two free parameters; Tfb, 
which sets the dissipative timescale, and n, which sets the effective equation of state 
in star-forming regions. Via this model, feedback energy can provide pressure sup¬ 
port to regions of gas that are thermally cold. Using a self-consistent disk galaxy, 
modeled after a local Sbc spiral, in both isolated and major-merger simulations, we 
investigate parameterizations of the feedback model that are selected with respect to 
the quiescent disk stability. These models produce a range of star formation histories 
for disks evolved in isolation, or duri ng a major merger, yet all are consistent with 
the star formation relation found by iKennicutO (|l998l) . We suggest that this result 
is produced by the adopted recipe for star formation and is not a byproduct of the 
feedback model. All major mergers produce a population of new stars that is highly 
centrally concentrated, demonstrating a distinct break in the surface density 
profile, consistent with previous findings. The half-mass radius and one-dimensional 
velocity dispersion are affected by the feedback model used. In tests with up to an 
order of magnitude higher resolution, the star formation history is nearly identical, 
suggesting that we have achieved a numerically converged star formation history. Fi¬ 
nally, we compare our results to those of previous simulations of star formation in 
disk-galaxy major mergers, addressing the effects of star- formation normalization, 
the version of smoothed particle hydrodynamics (SPH) employed, and assumptions 
about the interstellar medium. We conclude by suggesting several methods by which 
future studies may better constrain feedback models. 

Key words: galaxies: interactions - galaxies: evolution - galaxies: starburst - galax¬ 
ies: formation - methods: numerical. 


1 INTRODUCTION 

It is commonly assumed that galaxy interactions trig¬ 
ger powerful bursts of star formation. This assump¬ 
tion is supported by a variety of observational data 
which links morphologically disturbed or outright 
merging gal axies to enhance d levels of star forma¬ 
tion (e.g., Larson fc TinsleJ 


197^;_ Josegh_^^^ig^ 


19851: ike nnicutt et ahl Iiq'stI: iBarton Gillespie et al.l 20031 : 


Lambas et al. 1^0031 : Nikolic et al. 20041 ) . In fact, the most 
vigorously star-forming galaxies in the local universe, the 
ultra-lumino us infrared galaxies (U L IRGs), are nearl y 
all mergers (ISanders fc Mirabdl 1 19961 : iBorne et al.1 12OO0I) . 


Indeed, one of the fundamental assumptions of many galaxy- 
formation models is the efficient conversion of cold gas into 
new stars whenever two galaxies merge. These assumptions 
appear to be neces sary in order to reproduce the abun dance 
of Lyman-Break ( Somervil^,_Primack&^]abM 200l[) and 
submillimeter (iGuiderdoni et al. 19981 : Baugh et al. 120051) 
galaxies. 

The numerical simulations run to-date have supported 
the theory that disk-galaxy mergers yield a burst of star for¬ 
mation. The earliest work established that galaxy mergers 
can concentrate a large amount of gas in the central regions 
of the remnant dNegroponte fc Whitelll98J : lHernaui^ll989l : 
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iBarnes fc Hernauistlll99ll . Il996l ~). Later, by including sim¬ 
ple recipes to capt ure star formation a n d kinetic fee dback 
from massive stars, iMihos fc HernauistI (|l994cl . Il99fil . here¬ 
after MH96) produced the first quantitative results for the 
star formation induced during a major merger. These au¬ 
thors found that mergers between two equal mass gas-rich 
galaxies, so-called major mergers, trigger star formation at 
every close passage of the galactic nuclei, culminating in a 
burst at the final merger. In total, these starbursts convert 
~ 80% of the original gas mass into stars, an efficiency ap¬ 
parently independent of galaxy orientation or the presence 
of a spheroidal bulge in the primary disk galaxy. 

_ Subsequent work has become more ambitions. ISpringeil 

ll2000l. hereafter SOO) employed a slightly more elaborate 
treatment of the interstellar medium (ISM) by accurately 
treating gas heating and cooling. More recentl y, simulations 
have included shock-induced star formation IIBarn^ 2004 1 


and ac creting black holes liSpringel, Di Matteo fc HernauistI 
l2005bl '). While all of these studies have supported the fact 
that major mergers induce star formation, it has become 
clear that the star-formation efficiency is highly dependent 
upon the star formation and feedback assumptions. Because 
we lack a detailed theory of star formation and the influence 
young stars have on their surroundings, there are no clear 
means to include these physical processes in numerical sim¬ 
ulations. The typical approach is to formulate physically or 
observationally motivated recipes whose impact is controlled 
by one or more free parameters. This methodology is also 
motivated by the finite numerical resolution achievable. Cur¬ 
rent state-of-the-art numerical simulations evolve individual 
fluid elements which each represent ~ 10® Mq and 100 pc. 
As an example, star formation on kpc scales is observation- 
ally found to trace the gas density in what i s known as the 
Schmidt law dSchmidtll 19591 : iKemiicut^l 1998ll . Motivated by 
this relation, simulations can employ a recipe that converts 
gas to stars based on the local gas density. In this sense, 
a simulation can capture the sub-resolution star formation 
averaged over scales that are numerically resolved. Unfortu¬ 
nately, there exists no clear analog to the Schmidt-law for 
the effects of massive stellar winds and supernovae. 

Effective parameterizations for stellar feedback have 
been included in numerical simulations for the past decade. 
The first attempts used “thermal feedback,” in which en¬ 
ergy released by supernovae simply increases the thermal 
energy of the gas (iKat jll9^ 1 . Subsequent tests have yielded 
one of the most robust results regarding feedback: ther¬ 
mal feedback is ineffective at regulating star formation due 
to the short cooling times (and ra pid thermal energy loss) 
in dense star-forming regi ons (SOO; iNavarro fc Whitgll993l : 
iMihos fc Hernauisdll994m . Extensions of thermal feedback, 
in which a fraction of the super nova energy is imparted 
as kinetic energy to near by gas (INavarro fc White! 1 19931: 
IMihos fc Hernauisdll994^ , have proved effective at produc¬ 
ing stable, constantly star-forming disk galaxies, yet have 
the negative property that their implementations are reso¬ 
lution dependent (SOO) and they have difficulty stabilizing 
disks with a large gas fraction. Beyond this kinetic feed¬ 
back, most other recipes have loo ked to prevent t h e effi - 
cient radiation of feedback energy. iGerritsen fc Ick3 (Il997l l 
used UV radiation from young stars to offset efficient cool¬ 
ing, SOO created a turbulent reserv oir that retained feed- 
back energy for a set time-scale, and iThacker fc Couchm^ 


ll200iri simply inhibited cooling for 30 Myr subsequent to 
the in ject ion of feedback energy fsee lThacker fc Couchmar] 
I 2 OOOI and iKav et al.l[2002l for a comparison of various feed¬ 
back methods). More recently, methods have be en developed 
to model the multiphase nature of the ISM dYe nes et aP 
I 997 I: iHultman fc Pliaras^ll999l: ISemelin fc Combesll2002l : 


Marri fc White! I 2 OO 3 I : Smingel fc Hernauis'J 2003l j and the 

effects that supern ovae blastwaves may have on the ISM 
(IStinson et al.|[200^ 1 . While each of the above models has ad¬ 
vantages and disadvantages, many of these methods choose 
to regulate star formation in an equivalent manner: by us¬ 
ing feedback energy to effectively pressurize the ISM in star¬ 
forming regions. 

In this paper, we address the issue of supernova feed¬ 
back in more detail. Using a generalization of the star- 
formation and feedback model of SOO, we explore a range 
of parameters and determine how the resulting starburst is 
affected. This is a slightly different approach than previ¬ 
ous work, in that instead of physically motivating a specific 
feedback model, we attempt to understand how assumed 
feedback parameters affect star formation. Both the over¬ 
all efficiency of feedback and its dependence on gas density 
through an effective equation of state are explored. One of 
our long-term goals is to uniquely determine a single pa¬ 
rameter set, or feedback model, which will be used for a 
large series of major- and minor-merger simulations. Unfor¬ 
tunately, it is difficult to define a checklist against which 
each model can be measured, and hence it is difficult for us 
to endorse one model as significantly better than any other. 
We can, however, identify trends that are present in our var¬ 
ious models and use these to guide parameter exploration in 
future work. 

An important aspect of the work presented here will 
be to detail how assumptions made by MH96 and SOO af¬ 
fected their merger-induced starbursts. Even though these 
papers arrived at similar conclusions, the galaxy models, 
temperature and equation of state of the ISM, implementa¬ 
tion of feedback, and N-body/smoothed particle hydrody¬ 
namics (SPH) code were were all different. This last point 
is particularly relevant as the conservation of energy and 
entropy depends on the formulation of SPH as well as the 
numerical resolution ( Hernaui^ll993bl : iThacker et al.l[2000l : 


ISpringel fc HernauistI 2002l l. In general, all versions of SPH 
converge to an identical solution for high resolution, but 
can be quite diff erent at low resolution . This result has re¬ 
cently motivated ISpringel fc HernauistI ll2002|j to develop a 
new formulation of SPH which manifestly conserves both 
energy and entropy, when appropriate. We will investigate 
what effects this new version of SPH has and how this com¬ 
pares to older formulations of SPH. We also address the nu¬ 
merical convergence of our simulations by performing merg¬ 
ers with up to ten times the resolution. In additional tests, 
we vary the assumed gas metallicity and star-formation ef¬ 
ficiency and determine the resultant star formation. 

In our study of the star-formation response to differ¬ 
ent parameterizations of feedback, we hnd that the proper¬ 
ties of the merger remnants are also affected. This is sig¬ 
nificant because mergers between spiral galaxies are consid¬ 
ered the most likely mech anism to form galactic spheroids. 
This “merger hypothesis” (iToomre fc Toom relll972t ^omr^ 


1977|) i s generally supported b y simulations I Hernauislll992l . 
1993 d : [Naab fc Burkertll2003l . SOO), and here we show that 
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details of the feedback model will affect the ability to test 
the merger hypothesis. We reserve an exhaustive study of the 
merger remnants to future work, but the variation of feed¬ 
back parameters alters the basic properties of the merger 
remnant, and we document these effects here. 

In short, we address the following issues in this work: 

(i) How does the star-formation rate and the global gas 
consumption depend upon the feedback model and feedback 
parameters? 

(ii) Is the merger remnant affected by the choice of feed¬ 
back parameters? 

(iii) Are our models consistent with previous simulations 
of disk-galaxy interactions? Does using a modified version 
of SPH that manifestly conserves both energy and entropy, 
when appropriate, change the star-formation evolution? 


This paper is designed to be one in a series which com¬ 
piles the star-formation and remnant properties for a large 
number of galaxy merger simulations. Prior w ork has ana¬ 
lyzed the kinematics of the merger remnants llPekel et al.l 
l2005h . the influence of dust on the observe d properties of 
the interacting system ij.Ionsson et al.ll200fil~l. the shapes o f 
the stellar and dark matter remnants dNovak et al.| 2006| l. 
and t he dissipative origin of elliptical galaxies ( Dekel fc CoxI 
l200fil l. Future work will investigate the star formation in¬ 
duced during mergers between unequal mass galaxies (Cox 
et al., in preparation), propose a model that predicts the 
properties of merger remnants (Covington et ah, in prepa¬ 
ration), and compare the simulations to observed samples of 
interacting and merging galaxies using non-parametric mor¬ 
phological classihcation systems (Lotz et ah, in preparation). 

We organize this paper as follows. To begin, ^details 
the numerical techniques and the physical ingredients imple¬ 
mented into the simulations, and !j3]introduces a disk galaxy 
which will be used to investigate our feedback models. This 
disk galaxy is simulated in isolation in H3.2I and during a 
major merger in also shows the star-formation his¬ 

tory (3131) and takes a brief look at the properties of the 
merger remnant l il4.4ll . We discuss our results in ^ includ¬ 
ing a resolution study (3S3J, a comparison to previous work 
(353J, and the effects of varying the ISM metallicity 1 ^)5.311 . 
Finally, ^provides some concluding remarks and a prospec¬ 
tus for future work using the techniques we have developed 
in this paper. 


2 NUMERICAL SIMULATIONS 

All numerical simulations perf ormed in this work use the 
N-Bo dv/SPH code GADGET l|SDringel. Yoshida fc Whit3 
[looj). Specifically, we use a variant of the publically avail¬ 
able v ersion which has the “conser vative entropy” version of 
SPH dSoringel fc Hernauistl[2nn2h . and additional routines 
that track the radiative cooling of gas and star formation. 
This version of GADGET was obtained from Volker Springel 
in 2001. In addition, we implemented several new features 
ourselves; stellar feedback, metallicity-dependent cooling, 
and the ability of each gas particle to spawn multiple new 
stellar particles. The following subsections describe in more 
detail the physical processes included in our simulations. 


2.1 Cooling 

Due to radiative cooling, baryons fall into the cen¬ 
ters of dark-matter halos and achieve the cold, dense 
environment conducive to the formation of stars. The 
ra diative cooling rate Anet is comp uted as described 
in IKatz. Weinberg fc HernauistI l|l99fil l. where the gas is 
treated as a primordial plasma and the ionization states of 
H and He are explicitly tracked under the assumption of 
collisional ionization equilibrium. In practice, the radiative 
cooling rate Anet then simply a function of the local gas den¬ 
sity pgas and internal energy per unit mass u. However, the 
radiative cooling of astrophysical plasmas also highly depen¬ 
dent upon the gas metallicity, we also r an merger simulations 
using the tabulated cooling curves of [Sutherland fc Dopital 
(|l993h . which allow for a range of metallicities. The results 
of these simulations are presented in il5.3l In all simulations 
the spatial and temporal changes in metallicity are ignored, 
although we keep track of these for us e in radiative trans¬ 
fer calculation s through the dusty ISM (|,Ionsson et al ] l200fil: 
|jonssonll200dj . Lastly, we remind the reader that, without 
adiabatic processes, the radiative cooling described here, 
in which molecular cooling is neglected, effectively sets the 
minimum gas temperature to 10‘^K. 


2.2 Star Formation 


All numerical simulations presented here include star for¬ 
mation. Star formation is assumed to be proportional to the 
local gas density and inversely proportional to the local dy¬ 
namical time-scale 


dp* Pgas 

dt td 


where c* is a free parameter determining the efficiency 
of star formation, and tdyn=(47rGpgas)~^^^. This star- 
formation recipe mirrors the approach comm only imple¬ 
mented in numerical simu lations (MH96; SOO; lKa.tzlll992l : 
ISpringel fc HernauistI [ioO^ L however we will briefly revisit 
this in § 15.2.31 

The star-formation prescription set forth by Equa¬ 
tion © was origin ally motivated by star-forming regions 
in our own galaxy llSchmiddll959l l but has been shown to 
hold in a wide range o f environments. The observational 
work of lKennicn td (1 199811 has shown that, on kpc scales, the 
aperture-averaged star-formation rate per unit area is corre¬ 
lated to the gas surface density. This relation is empirically 
determined to be 


Esfr = (2.5 ±0.7) X lO"'* 


M0pc- 


1.4±0.15 


Me 

yr kpc^ ’ 


( 2 ) 


and will be referred to as the “Kennicutt law” throughout 
the remainder of this paper. While the index 1.4 resem¬ 
bles the p* ~ Pgas implied by Equation ([T}, we must admit 
that it is not immediately clear that the three-dimensional 
Schmidt-law, based upon the gas volume density, should be 
equivalent to the empirical two-dimensional Kennicutt-law, 
based upon the gas surface density. Nevertheless, the simu- 
lati ons do obey as shown in 3331 

iKenniciittl (Il998l l found the existence of a surface den¬ 
sity threshold at ~ 1 0 Me pc~^. below which star formation 
dramatically drops (iKennicuttI Il998l : [Martin fc KennicuttI 
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[20o3). To capture this effect in our numerical simulations, 
we introduce a volume-density threshold, pth- Gas which has 
a density larger than pth is eligible to form stars at a rate 
given by Equation 0. We adopt pth = 0.0171 Mq pc 
a value similar to that used in SOO (MH96 did not in¬ 
clude such a density threshold) and confirmed to reproduce 
the observational suppression of star formation below 10 
Mqpc”^. Interestingly, Equation © also appears to hold 
for azimuthally-averaged gas densities and star-formation 
rates in individual galaxies, although the scatter seems to 
be much larger and the slo pe shallower dMartin fc KennicuttI 
l200ll:IWong fc Blitj[200l) . 

Our implementation of star formation is fully deter¬ 
mined by the one free parameter in Equation 0, c*. To 
help us select a numerical value for c* we define the charac¬ 
teristic time-scale for star formation as 


U = 


Pga 


tdyn 

C* 


(3) 


Equation © relates the star-formation time-scale to the 
dynamical time by our one free parameter c*. Both time- 
scales become short at high densities. Observations indicate 
that the median gas consumption time-scale is much longer 
than the dynam ical time-scale, suggesting that c* <C 1 
llKennicutd ligga b After experimenting with various val¬ 
ues for c* (see iEij we found that 0.03 provided a good 
fit to Equation although we caution that the obser¬ 
vational scatter is quite large - almost an order of magni¬ 
tude higher or lower. Hence there could be a range of al¬ 
lowed values for c*, all consistent with the Kennicutt law. 
Our value of c* = 0.03 fixes the time-scale for star for¬ 
mation to be 1.1 gigayear (Gyr) at pth. Many previous 
studies also include star formation via Equation ([T}. These 


both hieher JXatJ 1992 

: Navarro & White! 19931: 

Kobavashil 

I 2 OO 4 I: iGovernato et al.l 

2004: Kawata & Gibson 

2 OO 5 II and 

lower (MH96; SOO; 

Sommer-Larsen. Gotz & Portinarill2003l: 


authors investigated a range of val ues for c* (iBekki fc Shioval 
1 19981 : iThacker fc CouchmarJ[2000h as we do in H5.2[ or use a 
hybri d model where Ct, can change with redshift or envi ron- 
ment dSommer-Larsen et akll^OOSl : lOkamoto et al.ll200^ . 

Very little is known about the requisite conditions en- 
abli ng star formation to occur. Previous simulation s (see 
e.g.. lKat3ll992l : lKravtsovll2003l : [^vernato et al.ll20()ll have 
included star formation criteria such as gas convergence, 
Jeans instability, and temperature thresholds. With the res¬ 
olution of current simulations, where individual particles 
represent ^ 10® Mq of gas, it is not clear that gas conver¬ 
gence should be required. The Jeans criterion is dependent 
on the numerical smoothing length, and hence has the un¬ 
desirable behavior that it is a function of numerical param¬ 
eters. Finally, as mentioned above, our treatment of cooling 
sets the minimum gas temperature at lO'^ K, thus we do 
not expect gas to approach physically reasonable tempera¬ 
tures for star formation (^ 50K). For these reasons we have 
not included any of the above criteria in the work presented 
here. 

As gas is converted into stars at the rate dictated by 
Equation 0, one must decide how to handle the newly 
formed stellar mass. Since it is not computationally feasi¬ 
ble to form individual solar-mass particles, numerical sim¬ 


ulations commonly treat SPH particles as “hybrid” parti- 
cles, representing both ga s and new stellar material (e.g. 
iMihos fc Hernauisd Il994blf . While this procedure is rela¬ 
tively simple to implement, it has the disadvantage that new 
stellar material is dynamically coupled to the gas. 

An alternative approach, suggested by several author s 
('e.g.. li^tJl992l : ISDringel fc Hernauistl20o4lBottemall2003l ~l. 
is to treat star formation as a stochastic process. Under this 
assumption each gas particle spawns a specified number Ng 
of equal-mass star particles. Then, if each gas particle forms 
a new star particle according to the probability 

P = Ag (l-e"*/**) , (4) 

the stochastic star-formation rate will equate to the designed 
continuous star-formation rate when averaged over the en¬ 
tire simulation time. We follow the above authors when im¬ 
plementing this procedure, and thus at each timestep ac¬ 
tively star-forming gas particles draw a random number be¬ 
tween zero and one. If this random number is less than P, 
a new collisionless particle, of mass m^as/Ng, is created. We 
note that rugas is the original gas particle mass and that the 
SPH particle mass decreases as individual star particles are 
generated. For the majority of the simulations used in this 
work, we select Ng = 2, although we have tried varying Ng 
up to 10 and found little difference in the results (see the 
end of E© for details). 

Throughout this work stars formed during the simula¬ 
tion will be called “new”, as opposed to “young”, stars. We 
make this choice because at any point during the simulation 
there exists a mixture of newly generated stars which have a 
distribution of stellar ages. While many of these particles are 
part of a young stellar population, a small percentage will 
have ages of several gigayears by the time the simulation is 
complete. 


2.3 Metal Enrichment 


Although it is not part of the analysis in the present pa¬ 
per, we have also included metal enrichment owing to stel¬ 
lar winds and supernovae ejecta. The enrichment, which as¬ 
sumes a yield of 0.02 per solar mass of stars formed and in¬ 
stantaneous recycling, is performed in a continuous fashion, 
based upon the instantaneous star-formation rate provided 
by Equation 0. The star-formation rate is calculated for 
each particle and metals are recycled within the same parti¬ 
cle, i.e., metals are not distributed among neighboring par¬ 
ticles nor is the diffusion of metals considered. These metals 
are carried by each gas particle and subsequent generations 
of new star particles adopt the metallicity of the parent gas 
particle at the time when they are born. Under these as¬ 
sumptions, the metallicity of individual gas particles can 
only increase with time. However, once the gas particle is 
fully transformed into stars, these metals are forever locked 
in the stellar component. We intend to implement a more 
realistic approach where star particles return a fraction of 
their mass and meta ls back to the inters t ellar medium in 
the future (see, e.g., [Mosconi et al.| |200ll: iTornatore et al.l 


|2004 lOkamoto et al'l l2005l : Scannanieco et al. 20061 . 2005l ~l . 
Finally, as noted in § 12.11 even though we track the metal¬ 
licity of the gaseous component, this information does not 
affect the cooling rate. 
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2.4 Feedback 


As discussed in the introduction, many of the methods by 
which feedback from young stars is incorporated into nu¬ 
merical simulations involve artificially pressurizing gas in 
high-density, star-forming regions. Within SPH, where par¬ 
ticles represent fluid elements, this pressurization is imple¬ 
mented on a particle by particle basis. This general approach 
is also used for the present study but, because there exists no 
clear theoretical understanding of how feedback pressuriza¬ 
tion should depend on the physical state of the star-forming 
regions, we instead adopt a general formalism and inves¬ 
tigate the star-formation response to a variety of different 
parameter choices. 

As the starting point for our exploration into feedback, 
we begin by adopting the feedback model put forth by SOO. 
The relative simplicity of this model and the ease with which 
it can be implemented into GADGET motivate this choice. 
Within this model, and adopting the nomenclature of SOO, 
gas is assumed to consist of two energy reservoirs. One, the 
standard thermal energy u, can radiatively cool at a rate 
given by Anet- The second, q, is an “ad-hoc” feedback reser¬ 
voir. The pressure of this medium is then calculated as 


Peff = ( 7 - 1)(m + '7)p, 


(5) 


where p is the gas density (the gas subscript is henceforth 
dropped), and 7 is 5/3, consistent with the assumption of 
a monatomic ideal gas. Equation © demonstrates that re¬ 
gions which are thermally cold, possibly because the den¬ 
sity is such that cooling is very efficient, can be pressure 
supported if the feedback reservoir u. 

Feedback energy is assumed to come from core col- 
lapse supernovae and is calculated in a manner similar to 
iNavarro fc Whi^ (Il993l l. New stars form with a power-law 
initial mass function of slope -1.5 between 0.1 and 40 Mq, 
and every star above 8 Mq goes supernova releasing 10 ®^ 
ergs of energy. The short lifetime of these high-mass stars 
is ignored and supernova energy is instantly returned to the 
system. We have investigated the effects of including realis¬ 
tic delays in supernova feedback, but it made little difference 
to the results reported here. We also investigated a limited 
range of excursions from the initial-mass function and high- 
mass cutoff listed above. In practice, though, these changes 
can be compensated with variations in our star formation 
or feedback free parameters and we do not consider this 
further. The supernova feedback energy per solar mass of 
newly formed stars esN is assumed to be first released into 
the feedback reservoir q, i.e. 


f —^ =esNi — 

\ dt / feedback P dt 


( 6 ) 


The feedback reservoir subsequently thermalizes its energy 
as 



thermalize 



thermalize 




(7) 


where f{p) parameterizes the efficiency of thermalization. 

The original work of SOO chose f{p) = Pl^/p, i.e., the 
thermalization of feedback energy is slower in high-density 
regions. For star-forming gas whose pressure is determined 
by the feedback reservoir, this choice results in a “stiff” 
equation of state P ~ SOO provided two motivating fac¬ 
tors for this model. First, a self-gravitating sheet of gas with 


this equation of state has a vertical scale height that is inde¬ 
pendent of surface density, similar to observations. Second, 
if the star-formation law is given by Equation ©, then this 
same self-gravitating sheet has a star formation rate per unit 
area Esfr that scales as the gas surface mass density a to 
the 1.5 power, i.e., Esfr oc ct^ ®. In other words, the star for¬ 
mation expected for an isolated disk manifestly scales in a 
similar fashion to the observed Kennicutt law (Equation [2]). 

While the formulation set forth by SOO was quite suc¬ 
cessful at regulating star formation in both isolated and 
merging galaxies and, as designed, reproduced the observed 
star formation laws, there are several reasons to explore al¬ 
ternative formulations. First, it is not clear that all observa¬ 
tions have converged on vertical scale heights of disk galaxies 
that are independent of radius. Second, the observed star- 
formation law holds in physical environments, such as the 
centers of spiral galaxies and in merging galaxies, that are 
very different from isolated disks. Thus, while the SOO equa¬ 
tion of state naturally predicts isolated disks that follow 
Kennicutt, this may be a necessary, but not a sufficient con¬ 
dition. Lastly, a quick survey of the literature returns a num¬ 
ber of alternate formations of feedback which do not lead to 
the SOO equation of state yet are still able to regulate star 
formation consistent with the Kennicutt relation (see e.g., 
ISpringel et al.ll2005bll . 

For these reasons, we depart from SOO at this point and 
assume /(p) is a general power-law function of density, i.e., 


/(P) 


1 

T-fb 



(l-n)/2 

, forp > Pth 

- 1/2 

, forp < Pth 


( 8 ) 


where pth is the critical density for star formation and we 
have introduced two free parameters: rfb, the timescale with 
which feedback energy is dissipated (and hence the efficiency 
of the feedback), and n which sets the polytropic index. The 
power-law below pth is fixed at n = 2 to ensure that the feed¬ 
back dissipation timescale gets shorter for low-density gas. 
We note that setting rfb = y/pfh/P and n = 2 transforms 
our model to the SOO model. 

In order to gain some insight into our model, as well as 
facilitate the selection of parameters to explore, we next cast 
our model in terms of an effective temperature and pressure 
in star-forming regions. Specifically, we consider a region of 
gas where the adiabatic energy changes are negligible, but 
which is also at a density greater than pth and hence is ac¬ 
tively forming stars. Under these conditions the feedback 
reservoir will reach a steady state in which the energy input 
due to supernovae equals the energy lost to thermaliza¬ 
tion 0 , which is subsequently radiated away, viz.. 


esN 


fdyn /(p) 


esNC* 



/ \ G-i)/2 


(9) 


This results is star-forming regions that have an effective 
temperature and pressure given by 


Tefl ~ 


(7 - 1)P 

k 


esNC* 




( 10 ) 


and 


Pe« — (7 — l)esNC* 


/ Tfb \ 

■p(n+l)/2- 

Vldyn J 

- Pth 


( 11 ) 
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Table 1. Star-formation and feedback parameter sets used in ei¬ 
ther isolated galaxy or major-merger simulations, or both, c* is 
the star-formation efficiency and n the density dependence of the 
feedback energy thermalization timescale, is the thermaliza- 
tion timescale, and Teff the effective temperature, in Kelvin, at 
the threshold density for star formation SOO used the same 
feedback parameters as our n = 2 “low-feedback” model, but 
a lower star-formation efficiency which is below the Kennicutt 
value. 


Model 

c* 

n 

rfb (Myr) 

Teff at Pth (K) 

SOO 

0.004 

2 

0.83 

3.5 X 10® 

low 

0.03 

0 ,1,2 

0.83 

2.6 X lO'* 

medium 

0.03 

0 ,1,2 

8.3 

2.6 X 10® 

high 

0.03 

0 ,1,2 

82.7 

2.6 X 10® 


where JI is the mean molecular weight, and k is Boltzmann’s 
constant. 

As mentioned previously, and shown explicitly by Equa¬ 
tion dni, n fixes the polytropic index, i.e., the power-law de¬ 
pendence of pressure on gas density P oc (remember 

that tdyn oc p“^Z2) When n = 0, the pressure increases lin¬ 
early with the density, as in an isothermal gas. When n > 0, 
the pressure increases faster than linear, and in this case the 
gas has a “stiff” equation of state. In this work we will in¬ 
vestigate three choices of n: 0, 1 and 2, thus from isothermal 
to the same n = 2 as investigated in SOO. While the present 
study assumes that n is a constant, there is no reason why 
n could not be a functio n of gas density . In fac t, this is the 
case in a recent study bv ISoringel et ah (l2005bli where their 
“multi-phase” model results in n = 2.8 near pth and then 
gradually softens to n = 0.7 at higher densities. 

To minimize the effects of the feedback reservoir q in 
low-density regions, i.e. below the critical threshold for star 
formation, /(p) is fixed to (he-i n = 2) for p < pth- 

The coefficient rfb ensures that /(p) is continuous at pth. 

The numerical implementation of these star formation 
and feedback processes is handled in a continuous fashion on 
a particle-by-particle basis. At each timestep, active gas par¬ 
ticles can form stars and instantaneously deposit energy into 
the ISM. Thus, Equations m and 0 become source terms 
in the standard hydrodynamic equations that are evolved 
within GADGET. 


2.5 Parameter Sets 

The star-formation and feedback model presented in the pre¬ 
vious sections has three free parameters: c*, controlling the 
efficiency of star formation, rfb, determining the timescale 
of feedback thermalization, and n, specifying the feedback 
density dependence, or, as Equation m demonstrates, the 
equation of state. While the star-formation efficiency c* was 
fixed in §[I1 Tfb and n are essentially free parameters. 

The approach of the work presented here is to survey 
a range of feedback parameters which sample, but not ex¬ 
haust, viable possibilities. Unfortunately, it is unclear how 
to judge the viability of feedback parameters. One possibil¬ 
ity that has been frequently used in the literature is to re¬ 
quire that the feedback be sufficient to stabilize the gaseous 
disk of a quiescent spiral galaxy against dynamic instabili¬ 


ties. This appears to be a zeroth-order requirement because 
observed spiral galaxies are dynamically stable. While this 
procedure is well defined, it does not generally yield a unique 
parameter set, but rather sets a minimum amount of feed¬ 
back necessary to stabilize a specific galaxy model. Galaxy 
models that have high gas fractions, high baryon fractions, 
and slowly rising rotation curves will be more susceptible to 
disk instabilities and thus require more feedback. 

In the current paper, we will also use the quiescent disk 
galaxy stability as a tool to investigate parameter choices. 
However, the model disk galaxy we use is designed to repre¬ 
sent a Sbc Hubble-type galaxy (this will be introduced in the 
following section), and thus is a large, gas-rich spiral that 
represents one of the most gravitationally unstable galaxies 
observed. Since Sbc spirals are stable, it seems reasonable to 
require that any plausible feedback model will need to sta¬ 
bilize our Sbc model. In practice, the feedback model may 
be more stable, and thus we will investigate a range of pa¬ 
rameter sets and attempt to understand their consequences. 

In order to define a reference feedback timescale rfb, we 
look at the stability of our isolated gala xy model gas disk as 
meas ured by the Toomre Q parameter llBinnev fc Tremaine! 
1 1981 eq. 6-49), i.e., 


Q 


CsK 

TfGSgas 


( 12 ) 


where c^='yPp~^ is the sound speed, k the epicyclic fre¬ 
quency, and Egas the gas surface density. For our model, 
the pressure of star-forming gas, and hence its sound speed, 
is approximated by Equation (HU). We also note that ob¬ 
servations suggest that most disk galaxies reside very clos e 
to Q ~1 (iMartin &: KennicuttI 1200 ll : I Wong fc Blit j [200^ 1. 
Thus, by inserting HU into HU, and setting Q = 1 we 
arrive at an expression for the feedback efficiency 


tdyn 7 Sgas7rG \ Y p 

7(7-l)esNC*V K J \pth J 


(13) 


Thus, via Equation dSJ, Tfb can be determined from the 
star-formation efficiency c*, the feedback density depen¬ 
dence n, the dynamical time at a specified density and a 
disk-galaxy model which specifies k and E. Equation US 
has the nice property that if it is evaluated at p = pth, as 
we will do here, then the dependence on n is removed. In 
practice, both k and E are functions of disk radius, but their 
ratio reaches a minimum at approximately twice the expo¬ 
nential disk scale length and we adopt this value for use in 
Equation (HI. 

We will label the Q = 1 case as the “medium” parame¬ 
ter set and use this as a reference for explorations to a super¬ 
stable “high” parameter set and an unstable “low” param¬ 
eter set, which are ten times higher and lower, respectively, 
and should maintain the disk at Q ~ 3 and Q ~ 0.3. These 
values of rfb serve to sample the range of plausible feedback 
efficiencies, and our model disk galaxy will be simulated in 
isolation and during a major merger for each parameter set. 
Table [T] fists the three values of ra, as well as the labels that 
will be used to reference each parameter set from this point 
forward. 
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Table 2. Properties of the Sbc disk-galaxy model used in this 
work. Mvir is the virial mass, including dark matter and baryons. 
The dark-matter halo has a NFW profile with a concentration c 
and spin parameter A, before this is contracted due to the presence 
of the baryons. The baryonic mass is distributed into one of three 
components: an exponential stellar disk, an exponential gas disk, 
or a stellar bulge. The mass fraction of each baryonic component 
is given by rrid, m^, and m^. is the stellar disk exponential 
scale length and zq is the vertical scale height. The gas disk has a 
scale length a times R^. Both gas and stellar disks are rotationally 
supported with angular momentum and given in terms of 
the total angular momentum. Finally, the bulge has a spherical 
exponential distribution with scale radius Rf,. 


Sbc 


Virial Mass, Mvir 

8.12 xl0“ Mq 

Concentration, c=Rvir/7’s 

11 

Spin Parameter, A 

0.05 

Disk Mass Fraction, md 

0.048 

Radial Disk Scale Length, R^ 

5.5 kpc 

Vertical Scale Height, zq 

1.0 kpc 

Stellar Disk Spin Fraction, 

0.058 

Gas Mass Fraction, 

0.065 

Gas Scale Length Multiplier, a 

3.0 

Gas Spin Fraction, 

0.219 

Bulge Mass Fraction, m;, 

0.012 

Bulge 3D Scale Length, Rf, 

0.45 kpc 

^halo 

100,000 

^disk 

30,000 

Hgas 

30,000 

^bulge 

10,000 


3 ISOLATED DISK GALAXIES 
3.1 The Basics of Model Disk Galaxies 


In this section we describe the basic ingredients of our model 
disk galaxy. In general , these models are similar to thos e 
described iniHernquistI (IQQSal) . iMihos fc Hernaui^ (|l996l) . 
ISpringel fc Whit^ 199£ b and SOO, as they contain a stel¬ 
lar and gaseous disk, which is rotationally supported, a 
spheroidal stellar bulge, all surrounded by a massive dark- 
matter halo. In the following subsection we will describe the 
parameter choices for our fiducial model, which is designed 
to approximate a large, gas-rich Sbc disk galaxy. In § 15.2.11 
we will also build an exact replica of the model used by SOO, 
which will aid in our comparison; however, we postpone dis¬ 
cussion of this model until then. 

The initial disk is rotationally supported and is assumed 
to have the following distribution 

Md 


Pd{R,z) = 


A'KzoR'i 


sech 


f z \ 1 


[fTo) ( 

^ Rd) 


(14) 


where Md is the disk mass, Rd is the radial disk scale length, 
and 20 is the vertical scale height. The disk is in general 
composed of both gas and stars, with a gas mass fraction 


Md,gas 


(15) 


We also allow for the possibility that the gas and stellar disks 
may have different radial distributions. While we continue 
to assume that both distributions follow Equation m , we 
will set the gas radial disk scale length Rd,gas = ctRd- 


We also define md and mg, the stellar and gaseous disk 
mass fractions, as 


md 


Md ,stars 

Mvir 


and 



(16) 


(17) 


where Mvir is the total virial mass of the galaxy, including 
all components, and Md = Md,stars + Md,gas. 

The vertical structure of the stellar disk is fixed by re¬ 
quiring that the disk be Toomre stable, i.e., Q ^ 1. In prac¬ 
tice, this results in a 2o that is set to 20% of the radial 
disk scale length Rd, and a vertical velocity dispersion of 
~ 4 0 km s~^at Rd , both consistent with observations (see 
e.g., lBottemalll993fl . 

The vertical structure of the gaseous disk is more com¬ 
plicated as the temperature is fixed by the equation of state, 
rather than the velocity dispersion (which is, in turn, fixed 
by our c hoice of 2o). This form alism has recently been up¬ 
dated bv ISpringel et all (l2005bll to account for pressure gra¬ 
dients in the gas disk and hence, they claim, can be used to 
construct purely gas disks. While we have not thoroughly 
investigated the effects of their new method, nor are we us¬ 
ing it here, we note that after an initial equilibration period 
our disk, with a relatively high gas fraction, is quite stable. 
This will be shown in the following section. 

We also allow for the inclusion of a centrally con¬ 
centrated stellar bulge. The bulge, when included, is non¬ 
rotating and has a spherically symmetric exponential distri¬ 
bution 



The exponential distribution is supported by observations 
which show that late-type spirals typically have bulge pro¬ 
files more consistent with exponential than the comm only 
assumed profile Jde Jong||l996l : lBalcells et'al]|2003h . In 
order to allow for this sec ond possibi l ity, w e also allow our 
model bulges to follow the iHernauistI (Il990ll distribution 


Pb{r) 


Ml, a 
2n r{r -\- a)^ 


(19) 


We found very little difference in the star formation between 
the two bulge distributions, and thus all galaxies in this work 
contain a bulge described by Equation (na. 

The disk-bulge system is surrounded by a massive 
dark-matter halo. The dark matter profile is assumed to 
follow the Navarro-Frenk-Wh ite (NFW) fitting formula 
llNavarro. Frenk fc Whit3ll997l l , characterized by two quan¬ 
tities, Mvir, the total mass, and the halo concentration 
Rvir/rs, where Rvir is the virial radius and Va is a scale ra¬ 
dius. In all of the galaxy models studies in this paper, the 
dark halo is con tracted due to the adiabatic growth of the 
central b aryon s jBlumenthal et al.| ll984l : iFlores et al.lll99^ : 
IMo et alii 19981 : iGnedin et al.l 2004li . 


The angular moment um of the dark h alo is set by the 
spin parameter A. In the IMo et aP lll998ll model for disk 
formation, the specific angular momentum of the baryonic 
component is conserved, and thus fixing jd, the disk angular 
momentum fraction, sets the size of the corresponding disk. 
We keep this formalism for the stellar disk, but treat the 
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gaseous disk separately. The gas disk, whose size is a specific 
fraction (a) of the stellar disk, has a much larger angular 
momentum fraction, jg. 

Once the above parameters are selected, the compound 
alaxy i s constructed in a manne r suggested by iHernauisl 
199381 1. ISpringel fc Whit3 (Il999l l. and SOO. In short, this 
method fixes the mass distribution according to the above fo- 
mula, and then assumes the velocities can be approximated 
by the Jeans equations. We refer the interested reader to 
these works for more details. 

Each disk galaxy is represented by N particles. In i)5.1l 
we investigate the choice of N by running some higher reso¬ 
lution simulations. Motivated by the criteria of I Power et al.l 
(l2003ll . all simulations presented here adopt a gravitational 
softening length h = 400 pc for the dark matter particles 
and 100 pc for the stellar and gas particles. We remind the 
reader that, in GADGET, forces between neighboring par¬ 
ticle become non-Newtonian for separations < 2.3 h. 

Now that we have outlined the basic ingredients of our 
disk galaxy model, we next introduce the specihc model used 
for our investigation of star formation and feedback. 

3.1.1 Sbc Galaxy Model 

Since our primary goal in this work is to investigate star for¬ 
mation and feedback parameters, we select a galaxy model 
which is gas-rich and thus likely to be gravitationally unsta¬ 
ble. To this end we choose a large gas-rich disk galaxy whose 
properties are motivated by local Sbc galaxies. This type of 
initial condition has the added benefit that the large gas 
supply provides a signihcant fuel supply and is thus likely 
to generate a large burst of star formation when participa¬ 
tion in a galaxy merger. Furthermore, large, gas-rich disks 
such as this may be similar to disks at high redshift where 
gas fractions are presumably higher. 

Because we choose the parameters of this model from 
observations of local Sbc’s, we will need to make some in¬ 
ferences to translate the observed quantities to our model 
inputs. To begin, we note that the optical radius Ropt, dy¬ 
namical mass Mdyn, ratio of HI to Mdyn within the optical 
radius, and absolute blue lumino sity Lb are neatly compile d 
as a function of Hubble type in [Roberts fc HavnesI (Il994h . 
Thus, for instance, the radial disk scale length is determined 
from the optical radius by assuming the relation Ropt = 3.2 
X Rd, which is valid for a Freeman disk. 

Next, the bulge-to-disk ratio and bulge size are con- 
strained to be consistent with the Sbc galaxies as found in 
Ide Jonj (Il996l b Because we have selected a late-type spiral, 
this necessarily results in a small compact bulge. We can 
determine the total stellar mass by using the average Sbc 
B-band luminosity Ls and the B-band mass-to-light ratio 
of 1.0 found for late-type spirals dBell fc de JondboOlll . The 
total stellar mass and the bulge-to-disk ratio can be used to 
determine the exact mass in both the stellar disk and bulge 
components. We note that with these assumptions the bulge 
makes up slightly less than 10% of the baryonic mass. 

Sbc galaxies are gas-rich. Observationally, there is a 
wide range of gas properties associated with Sbc galax¬ 
ies, but in general, the gas component is roughly exponen- 
tially distributed and more ex tended than the stellar disk 
llBroeils fc van Woerdenlll994l l. Taking this into considera¬ 
tion, we assume the gas disk to have an exponential prohle 


with a gas scale radius Rg which is a factor a times the stel¬ 
lar disk scale radius, Rd. Based on t he most extended distri¬ 
bution for Sbc galaxies reported in iBroeils fc van Woerd^ 
lll994h . we choose a = 3. The mass of the gaseous com¬ 
ponents can be found from the ratio of HI to Mdyn within 
the optical radius, with our assumed profiles, corrected for 
helium. 

At this point, the mass and distribution of the bary¬ 
onic components are completely specified. However, the dark 
matter distribution is not. In order to get some handle on 
the dark mass and profile we can use the dynamical mass. 
In practice, we assume an adiabatically contracted NEW 
halo with a concentration of 11 and then iteratively search 
for dark matter mass which gives us the correct dynamical 
mass at Ropt. This procedure results in a dark halo which is 

8.1 X 10*^^ Mq. The resulting rotation curve has a maximum 
speed of approximately 204 km s“^ which is larger than the 
~ 19 0 km s~^ suggested by the baryonic Tully-Fisher rela¬ 
tion jBell fc de Jong|l200ll l. The large rotation velocity of 
this model is a typical bypr oduct of adiabat ically contract¬ 
ing an NFW halo (see, e.g., iGole et ^l2000ll . 

The above procedure does not use the halo spin param¬ 
eter to set the size of the stellar disk. Arguably, the actual 
scale lengths of galaxies are more well constrained than halo 
spin parameters. This also means that the angular momen¬ 
tum of the disk components is independent of the halo an¬ 
gular momentum. The halo spin param eter was set to the 
median value from iBullock et ^ (l200lll . Decoupling of the 
halo and disk angular momenta is also motivated by the dis¬ 
crepancy between the specihc angular momentum distribu¬ 
tion of a typical dark mat ter halo and that of an exponential 
disk jBullock et ahllioolh . 

Each disk galaxy is represented by 170,000 N-body 
particles, 100,000 of which represent the dark-matter halo, 
30,000 each for the stellar and gas disks, and 10,000 for the 
bulge. In H5.ll we investigate our choice of particle number. 
Tabled lists the properties for our Sbc disk-galaxy model. 

3.2 Star Formation 

The feedback model outlined in jdwas designed to stabilize 
the gas disk and prevent excessive star formation due to in¬ 
stabilities in the gas disk. To test our model, the Sbc galaxy 
model was simulated in isolation for one Gyr, with the last 
nine parameter sets in Table [T] The resulting star-formation 
rates are shown in Fig. [T] Except for nOlow and nllow, all 
models have continuous star formation at a roughly constant 
rate. The spikes of large star formation in the hrst ~ 100 Myr 
of the simulation result from the initial model not being ex¬ 
actly in equilibrium. In particular, the feedback reservoir q is 
empty at the beginning of the simulation and thus provides 
no pressure support in the initial disk. The low-feedback pa¬ 
rameter sets have star-formation rates of over 10 MQyr”*^, 
the medium-feedback parameter sets 2-3 Mq yr“^, and the 
high-feedback parameter sets < 1 Mq yr“^. The equation of 
state index, n, has a relatively minor effect on the average 
star-formation rate of the isolated galaxies. For instance, run 
nOmed has an average star-formation rate of 3.3 MQyr“^, 
while the average of n2med is 3.0 Mq yr~^. 

Figure [5] displays the gas surface density for the mod¬ 
els. All models with constant star formation are qualita¬ 
tively similar. In the nOlow and nllow models, the gas disk 
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Figure 1. Star formation for an isolated Sbc disk galaxy. The left, middle and right show the runs for the n =0, 1 and 2, feedback 
models, respectively. The “low” feedback model has the highest star-formation rate and the “high” feedback model has the lowest rate. 

— 1/2 

For a given n, the star-formation rate scales roughly as 


has fragmented due to the growth of instabilities in the in¬ 
ner disk. The n2low model also looks fragmented, but its 
star formation rate is still regulated at a constant rate of 10 
Mq yr“^. This evolution is consistent with other studies that 
have evolved isolated disks with insufficient feedback fSOO ; 


Bottema 

20031: iRobertson et al.ll2004l: ISDrineel et alj|2005bl: 

Li et al.l 

20051. l2006tl. We also note that these instabilities 


were likely seede d by noise in the potential llToomre|ll98ll : 
lHernaui^ll993al l and because the initial disk is not in per¬ 
fect equilibrium. In the stable disks, this noise is then swing 
amplified to produce a spiral pattern. Representing the halo 
with a fixed potential or increasing the number of particles 
used for the halo both delayed the onset of gravitational 
instabilities or spiral structure. 

The gravitational instabilities are also manifested in the 
star-formation rates which are much larger than the stable 
disks’ and are not maintained at a constant value. This is 
especially apparent in the nOlow and nllow models in which 
run-away star formation progresses until a large fraction 
(~ 60%) of the original gas is consumed. In general, how¬ 
ever, relatively little of the original gas is converted to stars 
during the 1 Gyr evolution we follow; only ~ 5% for the high 
feedback models and ~ 20% for the medium feedback mod¬ 
els. One reason for these low gas consumptions is the very 
extended distribution of gas. Half of the original gas mass 
is at radii greater than 27.6 kpc, whereas star formation 
predominantly takes place at the galaxy center and along 
transient spiral structures, always at radii less than 25 kpc. 
A more subtle feature is the transient spike of star formation 
at T ~ 0.05 Gyr, present in all models, that results because 
the initial system in not in perfect equilibrium. One rea¬ 
son for this is that the initial feedback reservoir is initially 
empty and requires a small amount of star formation to fill 
it sufficiently enough to provide pressure support. We tested 
one simulation that started with a filled feedback reservoir 
and while the initial transient star-formation rate was sup¬ 
pressed, the evolution beyond T = 0.1 Gyr was unchanged. 

Within our model, the time-averaged star-formation 
rate in a quiescent spiral disk is directly set by rfb. Figure [T] 
demonstrates that each low, med, and high feedback model 
forms stars at a nearly constant rate, regardless of n. We 
can understand the origin of this scaling by considering that 


the gas disk is in hydrostatic equilibrium within the galaxy 
potential. Since the gas is a small fraction of the total mass, 
the potential and thus the gas pressure (if it is in hydrostatic 
equilibrium) are expected to be fixed by the galaxy model. 
In other words, the gas pressure is constant regardless of 
the specific feedback model. For star-forming gas, where the 
effective pressure is provided by Equation GH), this yields 
a simple relation between the gas density and the feedback 
timescale 


P 


-l/(l + n/2) 

oc 


( 20 ) 


Since the integrated star formation rate is 

5'Fi?ocp*Focp*/pocp°®, (21) 

where the third proportionality comes from Equation ©, it 
follows that 

SFRQCT~^^^'^+"f ( 22 ) 


In the simplest, n = 0, case the star-forming gas has an 
isothermal equation of state and Equation (1221) predicts the 
star formation to scale as SFR oct^ .In our range of 
models, rfb was varied by factors of ten so that we expect 
the star formation rates to vary by a factor of \/To ~ 3. This 
appears to match the simulated star-formation rates quite 
well, where the values for the high, med, and low models 
are ~ 1, 3, and 10. 

These star-formation scalings are complicated by the 
distinct density threshold for star formation as well as the 
self-gravity of the gas. In particular, the n = 2 model ap- 
pears to follow the ' scaling the closest, yet its expected 
scaling is a much softer . The reason for this discrep¬ 

ancy, as we will see more clearly in § 13.51 is that most of the 
gas within these models piles up at the threshold density for 
star formation pth and thus has the effective temperature 
at this density (Tes in Table©. Thus, the gas scales as an 
isothermal gas rather than n = 2. 


3.3 The Kennicutt Law 

Star formation in the simulations was des igned to be con- 
sistent with the empirical Kennicutt law (iKennicuttI Il998l . 
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Figure 2. Projected gas density for the Sbc disk galaxy model for each parameter set shown after 500 Myr, or ~ 3 rotation periods 
at the half-mass radius, of evolution. The “low” parameter sets have led to the break up of the gaseous disk and hence accelerated star 
formation as shown in Fig. [T] 


Equation [2}. The success of this is verified in Figure |3] where 
we show the star-formation rate and gas mass averaged in 
azimuthal annuli for the isolated Sbc galaxy simulations. 
The Kennicutt law is shown as a straight black line, which 
the simulations closely track. 

While Figure [3] demonstrates that all simulated disks 
track the Kennicutt law at high gas surface densities, the 
simulations that included a threshold density pth fall below 
the Kennicutt law at lower gas surface densitie s. As men¬ 
tioned previously, a similar feature is observed bv iKennicuttI 
lll998ll near ~ 10 kpc pc“^. Figure|3]also shows that higher 
feedback models fall below the Kennicutt law at higher gas 
surface densities. This is a result of the increased vertical 
scale height owing to the high feedback pressurization. 


To ensure that it is pth that sets the departure from the 
observed Kennicutt law we reran the three medium feedback 
models {nOmed,nlmed, and n2med) with pth = 0. Figure[3] 
shows that these simulations closely track the Kennicutt law 
across all gas densities. 

3.4 Disk Scale Heights 

The differences in star-formation rates between the low- 
, medium-, and high-feedback parameter sets is accompa¬ 
nied by corresponding differences in the scale height of the 
gas. The medium-feedback simulations have scale heights of 
~ 300 parsec at a radius of 5 kpc, increasing to 450 parsec at 
25 kpc where star formation ceases. The high-feedback simu- 
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Log s gas (M© pc"') 

Figure 3. Star-formation rate per unit area and gas surface den¬ 
sity averaged within azimuthal annuli for the fiducial isolated Sbc 
disk galaxy at T = 0.1 Gyr. Shown with different symbols are the 
low (asterisks), med (continuous line), and high (open triangles) 
parameter sets for each value of n (0 = blue, 1 = green, 2 = red). 
Also shown are three med feedback runs that did not include 
an explicit gas density threshold for star formation. The straight 
(black) line is the empirical Ken nicutt law. This p lot can be di¬ 
rectly compared to Figure 3 from iKennicutf] {lOOa). 


lations, in contrast, have scale heights of 450 parsec at 5 kpc 
radius, increasing to 500 parsec at 25 kpc. The n = 0 and 1 
low-feedback simulations do not have a well defined scale 
height because the unstable disk has fragmented. The n = 2 
low-feedback simulation has a scale height of ~ 200 parsec, 
near the numerical resolution of our simulations. 

For a given value of ru,, the vertical structure is in¬ 
dependent of the value of n. This fact, which may be sur¬ 
prising at hrst, arises because for all these models the feed¬ 
back is efficient enough to limit the maximum gas density 
to slightly above the threshold density while most of the gas 
is at densities slightly below the threshold, even in the star¬ 
forming regions of the disk. Since n describes the scaling of 
effective pressure with density above the threshold, and the 
star-forming gas is confined to a narrow range of densities 
above the threshold, it follows that n has only a very small 
influence on the structure of the isolated galaxies. (As men¬ 
tioned earlier, for gas below the threshold density a value of 
n = 2 is used in all cases.) 

The fact that the star formation comes from gas at the 
threshold density also explains how a ~ 50% difference in 
scale heights between the med- and high-feedback sets can 
explain a difference in star formation rate of more than a fac¬ 
tor of 3. Without a density threshold. Equation © predicts 
that the integrated star-formation rate should scale as , 
but at the threshold density the effective star-formation law 
is much steeper than indicated by Equation ©. 


3.5 Phase Diagram 

To illustrate the effects of our pressurizing feedback, and 
its dependence on the free parameters rfb and n, we present 
temperature-gas density phase diagrams in FigureH] On the 
vertical axis of the phase diagram is the effective tempera¬ 
ture, which includes a contribution, if any, from the feedback 
reservoir. There are several lines plotted in Figure|T]that de¬ 
lineate distinct phases of gas. The horizontal line at 10“^ K 
shows the initial gas temperature as well as the minimum of 
the cooling function. Most low-density gas resides near this 
temperature. 

The feedback reservoir becomes signihcant for gas which 
is above the threshold density for star formation pth, shown 
by the vertical solid line in Figure [T] At high gas densities 
(Pgas ^ Pth) star formation acts to maintain the gas at a 
temperature and pressure specihed by Equations (nsi) and 
dm, which is shown by a dotted line in Figure [T] From these 
relations it is apparent that the slope of the dotted line is 
specihed by n, and the normalization by the combination 
c*Tfb. For example, when n = 0, star-forming gas is isother¬ 
mal with an effective temperature of 2.5x 10'^ K, 2.5x 10® K, 
and 2.5x10® K (see Table 0 for the low, med, and high 
models, respectively. 

Figure 0] demonstrates the effects of the equation of 
state in star-forming regions set by n. For instance, the 
nOlow model is nearly identical to no feedback at all. Gas is 
always close to 10“^ K, and the disk is insufficiently pressur¬ 
ized and subsequently fragments (see Figure [2]). In the case 
of n2low, however, the steep equation of state above pth 
means that high-density gas eventually becomes pressurized 
and limits the amount of fragmentation. 

In general, the gas phases below the threshold density 
pth are very similar for each high, med, and low model, 
regardless of n. For the medium- and high-feedback models, 
the high temperature and significant pressurization cause 
a large amount of the gas to pile up at pth. This effect is a 
byproduct of the long timescale for energy dissipation, set by 
Tfb, in the high and med models. In this scenario, because of 
recent star-formation, the feedback reservoir has been Hlled 
and the increased pressure support has pushed the particles 
to densities below pth. 

For the quiescent spiral galaxy we simulate here, we 
have not included any diffuse hot gaseous halo. Thus, 
gas only becomes hot (> 10^ K) because of star forma¬ 
tion. We also ex pect hot diffuse gas to be efficiently pro- 
duced by shocks (ICox et al.ll2004l ) and active galactic nuclei 
dCox et al.ll200d) that may attend the galaxy merger. This 
will be explicitly shown in S 14.31 


4 MAJOR MERGERS 

To generate a galaxy merger, identical copies of the iso¬ 
lated galaxy model are initialized on a prograde parabolic 
orbit with pericentric distance Rperi = 11 kpc (~ 2 Kd or 
~ 0.05 Rvir) and an initial separation of 200 kpc. One disk 
is in the orbital plane while the other is inclined by 30°. 
This generates a fast, nearly radial, collision consistent with 
orbit s found for dark-matter ha los in a cosmological simula¬ 
tion dKhochfar Burkerdl200^ ). 

Figure [5] displays the stellar (disk plus bulge) density 
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Figure 4. Each panel shows the effective temperature - gas density phase diagram for the isolated Sbc disk at T = 0.25 Gyr for the 
entire set of feedback models considered. The effective temperature is calculated as Tefr=^Peff/fcp, where Peff is given by Equation 
The vertical line shows the threshold gas density for star formation pthi S'lid the horizontal line at 1.2x10'* K is a guide to the eye that 
demarcates the floor of the cooling function. The expected effective temperature determined by our feedback model, Equation lllOll . is 
shown by a dotted line at densities above pth- 


projected onto the orbital plane and Figure [6] shows the 
projected gas density from an equivalent viewing angle. The 
simulation begins with a wide separation between the two 
galaxies, and they are essentially unaffected by each other’s 
presence. 

As the galaxies reach pericenter (T ~ 0.6 Gyr) each 
disk becomes tidally distorted. Long tidal tails carry out 
loosely bound material while the central regions form tran¬ 
sient bar-like structures (see T = 0.7—0.9 Gyr). At this point 
dynamical friction has increased the spin of the dark-matter 
halos and extracted orbital energy from the originally un¬ 
bound orbit. As the galaxies separate, a bridge forms that 
connects the two galactic centers, and at T = 1.1 Gyr the 
galaxies reach apocenter, with a separation of 80 kpc. 

The final merger is a messy process starting at T = 
1.6 Gyr and lasting for 200 Myr. The disks are effectively 
destroyed and we are left with a spheroidal-looking ob¬ 


ject (MH96; SOO: iHernauistI 0*99211 . Owing to gas shocking, 
kinetic energy is converted to thermal energy and subse¬ 
quently radiated away. This results in an offset between the 
stellar and gaseous components that torques the gas into the 
dense central regions of each galaxy (see T = 0.7 — 0.9 Gyr). 
Since we have tied the formation of stars to the gas den¬ 
sity thr ough Equation (llTl. a burst of star formation ensues 
(MH96; iMihos fc Hern(mistll994cl L FromT = 1.2—1.6 Gyr, 
the extended distribution of the gaseous component pro¬ 
duces long tidal tails while the stellar tails have already 
fallen back to the central galaxy. These gaseous tails are 
sharper and thinner than their stellar companions, and con¬ 
tinually rain back on the central disks. 

The final merger is particularly violent for the collisiona l 
gas as large amounts of gas shocking occur dCox et al.ll2004l i . 
These shocks provide a significant source of heating to the 
gas, and its evolution depends on its density. High-density 
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Figure 5. Projected stellar (disk plus bulge) density during the n2med major merger. Each panel measures 160 kpc on a side, and the 
number in the upper left of each frame is the simulation time measured in Gyr. The disk galaxy initially on the lower left is coplanar with 
the orbital plane while the upper right disk is tilted by 30°. The logarithmic density scale, in units of Mq pc“^, is given by a colorbar 
in the upper-left plot. 
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Figure 6 . Same as Fig. [5] except here the projected gas density is shown. The density scale is identical to Fig. [5] and is given in the 
top-left plot. 
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gas is able to radiatively cool on a short time-scale, and 
provides fuel to continue the starburst. Low-density gas, on 
the other hand, expands from the central region and con¬ 
tributes to a hot gaseous halo in the merger remnant. At 
T = 3.0 Gyr, 1.2 Gyr after the merger is complete, ~ 18% of 
the current gas content has cooled and has fo rmed a nuclear 
(< 10 kpc) disk, similar to what was found bv lBarn^ (l2002l l 
and lSpringel fc Hernaui^ (l2005li . The merger remnant also 
contains several tidal condensations, formed from loosely 
bound tid al material and composed entirely of gas an d new 
stars (c.f., iBarnes fc HernauisdllQQOl : IDuc et al.ll2004l ~l . 


4.1 Star Formation 

Figure [7] displays the star-formation rate during the sim¬ 
ulated major merger for all feedback parameter sets. The 
star-formation history is punctuated by several significant 
dynamical events. Prior to the first close passage {t < 0.5) 
the star-formation rate is that of two quiescent disks. The 
first encounter of the disks, at T = 0.58 Gyr, produces a 
strong shock between the gas disks and a ridge of high den¬ 
sity gas stretches across both disks. This feature generates 
an enhancement of the global star-formation rate, especially 
in the low and medium feedback models. 

Subsequent to the first passage, starting at 
T ~ 0.65 Gyr, tidally ejected gas is steadily stripped of 
its angular momentum and continuously falls back on the 
original disk galaxy. This process lasts for ~ 1 Gyr, but 
the resulting star-formation rate is strongly influenced by 
the feedback model. The low-feedback models provide gas 
with minimal pressure support. Gas quickly reaches high 
densities in the galaxy center, fueling a prolonged burst 
of star formation. High-feedback models provide much 
stronger pressure support, which limits the amount of gas 
at high density and hence the star-formation rate. 

At first glance, the prolonged star formation during and 
immediately after the first disk encounter seems at odds with 
the results of MH96 and SOO, whose progenitor disks, which 
contained a bulge as ours do here, demonstrated very minor 
star formation after the first encounter. In fact, the galaxy 
models employed by MH96 and SOO contained bulges more 
massive (one-third the disk mass) than that used here (one- 
tenth the disk mass), and a much smaller fraction of gas. 
As noted by MH96, this suggests that the suppression of 
gaseous inflow due to the presence of a bulge is very depen¬ 
dent upon the bulge mass and the disk composition. More 
work will be needed to quantify the relationship between 
bulge mass, galaxy composition, and early gas inflow lead¬ 
ing to bursts of star formation. 

The final merger is preceded by an increase in the star 
formation rate at T ~ 1.6 Gyr. At this time, gas begins to 
see the nuclei as a single potential, while it is also being 
stripped of angular momentum due to the abundant shocks 
that occur. At T = 1.8 Gyr, the final coalescence occurs. 
The strongest burst of star formation precedes the merger 
by ~ 50 Myr in every model. The strength and duration 
of this burst is a strong function of the feedback model. 
Models with isothermal equations of state have a strong, 
prolonged burst (except the nOlow and nllow models, which 
by this time have consumed a large portion of the available 
gas), while the models with stiffer equations of state suppress 


Table 3. Star formation properties of major merger simulations. 
Here e is the fraction of the original gas consumed by star forma¬ 
tion during the entire 3 Gyr simulation. SFRmax is the maximum 
star formation rate during the merger. Tlirg is the duration, in 
Gyr, at which the star formation rate is greater than 20 Mq yr“t^ 
corresponding to when this galaxy pair would be classied as a lu¬ 
minous infrared galaxy, LIRG. 


Model 

e 

X 7 

nj 1 

Tlirg 

(Gyr) 

nOlow 

0.70 

78 

1.55 

nOmed 

0.55 

159 

1.31 

nOhigh 

0.25 

51 

0.25 

nllow 

0.69 

84 

1.72 

nlmed 

0.51 

105 

1.29 

nlhigh 

0.24 

43 

0.08 

n2low 

0.67 

75 

1.85 

n2med 

0.47 

56 

1.34 

n2high 

0.18 

44 

0.04 


large bursts of star formation due to the efficient pressure 
support provided by the feedback reservoir. 

To provide a more quantitative view of the star forma¬ 
tion induced by the galaxy merger, we calculate the fraction 
e of the original gas converted to stars during the simulation. 
Table |3] lists e for each merger. Two trends are noticeable 
for e, both consistent with the isolated disk galaxy simula¬ 
tions: first, high-feedback models consume less than half as 
much of the original gas as the low-feedback models do, for 
equivalent n. Second, models with a stiffer equation of state 
(i.e., higher n) suppress gas consumption, although this is a 
much smaller effect. 

Another relevant measure of the global star formation is 
the peak star-formation rate SFRmax, which is also listed in 
Table |3] The equation of state n has a much larger influence 
on SFRmax than on e. For example nOmed has a SFRmax 
roughly three times larger than n2med, yet the overall gas 
consumption is larger by only 8%. Similarly, the peak star- 
formation rate of n2high is only 20% lower than n2med, 
while the gas consumption is less than half as much. 

The different dependence of e and SFRmax on the feed¬ 
back parameters results from the varying effects of Tfb and 
n. The effective pressure is proportional to Tfb (see Equa¬ 
tion [13, thus large values produce significant pressure sup¬ 
port and restrict the amount of gas that can reach star¬ 
forming densities. Thus the overall gas consumption e is 
primarily a result the value of Tfb. On the other hand, once 
gas reaches star-forming densities, n determines the abil¬ 
ity of this gas to reach even higher densities, and because 
the star-formation rate scales non-linearly with density (see 
Equation [TJ the maximum star-formation rate SFRmax is a 
strong function of n. 

It is interesting to note that the simple star-formation 
rate scaling, derived as a function of Tfb in Equation (I22II . 
that only moderately described the behaviour of the iso¬ 
lated galaxies, performs much better for the merging galax¬ 
ies star-formation rates. Recall that the expected scaling 
was SFR oc Applying this to the maximum star- 

formation rate, we would expect a ratio of 3.2, 2.2, and 1.8 
for the n = 0, 1, and 2, respectively. The ratios between the 
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Figure 7. Star formation for a major merger between two identical Sbc disk galaxies. The n = 0,1, and 2 feedback models are shown 
in the top, middle and bottom panels, respectively. In every instance the low feedback model has a higher star-formation rate than the 
high feedback model between T < 2 Gyr. The final merger occurs at T = 1.8 Gyr and this is coincident with the highest peak of star 
formation in all models with sufficient feedback to stabilize the isolated disks, i.e. the medium- and high-feedback parameter sets. 


med and high runs from Tableware 3.1, 2.4, and 1.3. Be¬ 
cause of the unregulated star formation in the low models, 
the scaling appears to break down for these. 

One interest in studying star formation induced by 
galaxy mergers is their possible relation to lum inous in¬ 
frared galaxies (LIRGs; Isanders fc Mirabell Il996l l. LIRGs 
are classified as galaxies which have an infrared luminos¬ 
ity greater than 10 ^^ Lq and hence are formi ng stars at 
> 20 Mq yr“^(using the L/_r-SFR conversion of iKennicutu 
llQQSt l. Listed in Table El is Tlirg, the total time during 
the merger in which this system would be detectable as a 
LIRG. The choice of feedback has a significant effect on the 
length of time this merger would classify as a LIRG. We 
note that we have not distinguished between the two pro¬ 
genitors, and thus do not account for viewing angles where 
the two galaxies would be classified as separate. Nor have 
we included a full treatment of the obscuring effects of dust 
and its re-radiation in the infrare d (but see Ijonsso nl l2004 
Ijonsson et al. I l2006l: [j onsso illooi). 

It is clear that these three quantities do not describe 
the full star-formation history of each model. Instead, they 
are gross measures, useful to characterize the effects n and 
Tfb have on the star formation induced by a galaxy major 
merger. Future work will compile quantities such as these for 
a wide range of initial conditions, assumed merger orbits, 
and merger mass ratios. These relations will be useful in 
analytic or semi-analytic models of galaxy formation, and 
in using observations to clarify which feedback models are 
realistic. 

As a final remark related to the star-forming histories 
shown in Figure [7] we draw attention to the oscillatory star 
formation that occurs in nearly every model, but is most 
pronounced in the high-feedback models. These oscillations 


occur primarily after peaks in the star formation rate and 
are thus clearly visible after the first passage and the final 
burst. However, even the isolated disk high-feedback mod¬ 
els contains some oscillatory star formation. The periods 
of oscillation are almost exactly 2 x Tfb, the thermalization 
timescale, indicating that the oscillations result from the in¬ 
terplay between pressurizing feedback efficiently quenching 
star formation followed by the dissipation of this pressure 
support, gravitational collapse, and an increase in the star- 
formation rate. 


4.2 The Kennicutt Law 

In § 13.31 we demonstrated the the individ ual spiral galax- 
ies matched the empirial Kennicutt law (iKennicuttI llQQSl , 
Equation 0, as they were designed to do. In Figure [S] we 
show the star-formation rate and gas mass averaged in a 
2 kpc aperture plotted at representative times during each 
merger. As was done previously, we overplot the Kennicutt 
law with a solid black line, and the observed envelope for 
which all observations are contained is plotted with a dotted 
line. Similar to the isolated galaxies, the merger simulations 
closely track the empirical star-formation law. There are a 
few points that are above the upper envelope, however these 
are from the unstable nOlow and nllow runs. The fact that 
all of the stable models agree with the Kennicutt law, re¬ 
gardless of n and rfb, suggests that it is our formulation of 
star formation according to Equation © and the value of 
the free parameter c* (which are identical for all runs) that 
determines the agreement with the Kennicutt law. 
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Figure 9. The effective phase diagram for the galaxy merger at T = 1.75 Gyr for all models. The calculation of Tefj and all of the lines 
are identical to Figure|4] The top-left panel {nOlow) denotes the three phases of gas, hot, cold, and star-forming as described in the text. 


4.3 Phase Diagram 

At this point we revisit the phase diagrams originally in¬ 
troduced in §[33] in the context of isolated Sbc disk galax¬ 
ies. We expect two differences when we look at the phase 
diagrams of the interacting galaxies in comparison to the 
isolated disks. First, the collisional shocks that attend the 
interaction generate a significant amount of hot low-density 
gas. Second, because of the interaction, much more gas will 
reach high densities. 

Figure [5] shows the phase diagram for the galaxy merger 
of each parameter set at a point in time 50 Myr prior to the 
final coalescence. At this point, most models are near their 
peak of star-formation. All gas resides in one of three distinct 
regions in the phase diagram: hot, cold, or star-forming. 

Gas which is below the critical density threshold for star 
formation pth is designated as either in a hot or cold phase, 
depending on whether its temperatnre is above or below 
1.2 X 10^ K, denoted by a horizontal line in FignrejOj When 
the simulation begins, as well as in the isolated disks, all non¬ 


star forming gas is cold owing to efficient cooling. This gas 
remains cold until the first passage, when gas contained in 
the tidal tails is ejected and adiabatically cools to ~ 10® K. 
Gas that is not part of the tidal material collides with gas 
from the other disk. These collisions generate shocks that 
heat the gas to temperatures well above 10“^ K where it 
is deemed hot. Because the cooling time of this hot, low- 
density gas is greater than the merger time-scale, most of 
this gas expands out of the central region and it remains in 
this phase until the merger is complete. 

During the final merger, the significant collis ion results 
in sh ock heating, and therefore a lot of hot gas iIGox et al.l 
|2004| . While most of the gas remains hot after the final co¬ 
alescence (T > 2.0), some of the shock-heated gas cools and 
flows to the galaxy center. This gas, as well as the contin¬ 
ual infall of tidal material, fuels the steadily decreasing star 
formation. 

Gas which is at, or above, the threshold density for star 
formation is strongly affected by the feedback model. In this 
sense, the merger is similar to the isolated disk. The main 
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Figure 8. Star-formation rate per unit area and gas surface den¬ 
sity averaged within an azimuthal aperture of radius 2.0 kpc, one 
point type for each simulation. For each run points are plotted 
every 200 Myr during the interaction, and every 50 Myr during 
the final merger. The solid line is the empirical Kennicutt law, 
and the dotted lines represent the envelope in which all observa- 
tion reside . This plot can be directly compared to Figure 6 from 
iKennicutJ |[l99^. 


difference is that gas in the merger reaches much high den¬ 
sities than the in the isolated case. Figure [U] shows that, like 
the isolated disks, the effective temperature of star-forming 
gas is close to Equation (nnj, which is delineated with a 
dotted line in the figure. 

Figure |9] also shows a large spread in the effective tem¬ 
perature of star-forming gas. This spread is considerably 
larger than in the isolated disks and is a direct result of the 
varied history of shock-heating, star formation, and feed¬ 
back that attends the violent nature of the merger. Figure[9] 
also demonstrates that the med- and high-feedback models 
have a significant amount of gas that piles up at pth owing 
to the pressure support provided by these models. This fea¬ 
ture also relates to the oscillatory star formation discussed 
at the end of the previous section as gas continually cycles 
into and out of star-forming regions based upon their recent 
star-formation history, current dynamics, and dissipation of 
feedback reservoirs q. 


4.4 Merger Remnants 

Most of the merger simulations run to-date have focused 
on the role disk-galaxy major mergers may play in forming 
ellipti cal g alaxies, as first suggested by iToomre fc Toornra 
and iToomrel dlQTif ). and commonly referred to as the 
“merger hypothesis”. To this end, studies have addressed 
whether galaxy mergers can reproduce the surface 

brigh tness profile which many elliptical galaxies appear to 


have |Barnesl 19921: Herno 

uisd 19921. 1993d: Hernauist et al.l 

19931: fMihos & Hernauisd 

1994al: [Barnes & Hernauisd 199(31: 

Soringell 2000l: Naab & Truiillq 2005|). the kinematic 


structure of merger remnants (iNaab fc BurkertI 200d; 
iGonzalez-Garcfa fc van Albadal l2005a T Gox et al. l2006l l. 
or if the merger remnants fall on the fundamental plane 
dBekki fc Shiova 1997 : Bekkil 199^: iNaab fc BurkerO 20031: 


Dantas_et_^ 20031: iNipoti et ^ l2003l i^ Bovlan-Kolchin et alJ 

20051. l2006l: iRobertson et al.ll2006ll . a tight correlation be¬ 


tween the velocity dispersion, effe ctive radius, and luminos¬ 
ity, on which a l l ellip ti cals reside |lDjorgoyski__&M)ayij 19^ ; 
Dressier et al.l Il987l: iBendy et al.l Il992l : iBernardi et all 


20031 : Padmanabhan et al. 2004). However, very few of 


these simulations have included recipes for star formation 
and feedback or explored how these parameters affect the 
remnant. In light of this, this section looks at the proper¬ 
ties of the merger remnants to determine if studies of this 
sort might be sensitive to the prescriptions used to include 
star formation and feedback. We reserve a complete anal¬ 
ysis of merger remnants for future work because this will 
require a much larger set of runs spanning a range of progen¬ 
itor galaxy masses. Here, we simply present the stellar mass 
profile, the remnant size, as measured by Re, the half-mass 
radius, and the one-dimensional central velocity dispersion 
for each merger remnant. We show that the size and velocity 
dispersion depend strongly on the feedback model while the 
stellar mass profile is only weakly dependent. 


4-4-1 Stellar Mass Profiles 

Figure [TO] shows the azimuthally averaged surface den¬ 
sity profile versus the fourth root of the projected radius 
for all merger simulations. A pure r^^*-type profile would 
be a straight line in this plot. Similar to the work of 
iMihos fc Hernauisd lll994ali and SOO, the merger remnants 
show central cusps that are disjoint from the outer density 
profiles. The central-density excess begins at R ~ 200 pc, 
near our resolution limit, and is comprised of both stars 
that are newly formed during the merger (shown as a dot¬ 
ted line), as well as the pre-existing bulge. In higher reso¬ 
lution simulations (see i)5.1ll the profile changed very little, 
as indicated by the overplotted profile in the n2med panel, 
suggesting that it is not just poor resolution which generates 
the central profile cusp. Further, because all models shown 
in Figure [TO] display a cuspy profile, it does not appear to be 
a result of the strength or density-dependence of the feed¬ 
back model, although there is a trend for the high feedback 
models to form fewer stars and thus the surface density is 
predominantly bulge-dominated toward the center. 

The ubiquitous nature of the central cusp indicates that 
it is a byproduct of our initial galaxy model. As noted by 
SOO, the remnant has not forgotten about the dense stel¬ 
lar bulge in the progenitor disk. Even a collisionless run, 
shown as a dashed line in the lower-right plot, demonstrates 
a break in the surface density. We also performed one simu¬ 
lation without a central bulge which produced a profile that 
exhibited a more subtle upturn toward the center. Not sur¬ 
prisingly, the central density excess in this case was solely 
composed of stars formed during in the simulation. T hese 
results are consistent with iMihos fc Hernauisj lll994a^ and 
SOO. 

It is not yet clear if the surface densities plotted in 
Figure [TO] are consistent with observations, or even if there 
exists a definitive sample with which to compare these pro¬ 
files. In a recent study of the K-band surface-brightness pro- 
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Figure 10. Surface density of each merger remnant projected onto the orbital plane versus the fourth root of the radius. Other projections 
are qualitatively similar. The solid line is for all stellar mass and the thin dotted line is stellar material formed during the major-merger 
simulation. The arrow in each plot points to the gravitational softening length of 0.1 kpc. The n2med plot, in the lowest middle panel, 
also shows the profile of an identical simulation with ten times the resolution and another where gas particles can spawn 10, instead of 
the standard 2, stellar particles. These variant simulations produce nearly identical profiles. The n2high plot, in the lower right-hand 
panel, contains a dashed line showing the total stellar profile of a purely collisionless run. 


files of 51 merger remnants bv iRothberg fc JosephI ll2004ll . 
16 (30%) had profiles steepe r than Si milar studies of 
ellipt icals llBvun et al.l Il996l i and ULIRGs dVeilleux et al.l 
l2002l i have found a much smaller percentage that deviate 
from indicating that the discrepancy between the ob¬ 
servations and simulations is far from obvious. In fact, dur¬ 
ing the one Gyr evolution subsequent to the final galaxy 
merger, the central excess slowly increases due to the rem¬ 
nant gas slowly being converted to stars. This trend might 
have been abated, or even reversed, if a more violent in¬ 
jection of energy would have dispersed, or removed, the re¬ 


maining gas thereby causing the density profile to remain 
constant or slightly decrease. While no source of energy is in¬ 
herent to our feedback model, recent work bv ISpriuTOl et al.l 
(l2005ah has shown that an active galactic nucleus (AGN) 
can indeed expel large amounts of gas and completely shut 
off star formation. Gorrelations between m ergers and AGN 
dVeilleux et_ajj []T9^ WvMie&Jjoeyi^OO^ and r nergem and 
outflows I Mart inl ll999l. l2005l: iRupke et ^l2005l. [20^ pro¬ 
vide evidence for this scenario. Moreover, merging black 
holes at the galactic center could also influence the cen- 
tral stellar distribution (see, e.g., iMilosavlievic fc Merritd 
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Table 4. Projected remnant properties averaged over 1000 ran¬ 
dom lines of sight. Re is the half-mass radius of all stars, while 
Re,os is the half-mass radius of the old stellar component, both 
disk and bulge, and Re,ns is the half-mass radius of stars gener¬ 
ated during the course of the merger simulation, a is the velocity 
dispersion of all the stars measured within half Re, and 0-5 is 
the velocity dispersion measured inside an aperture of 5 kpc. The 
total and new star effective radii are not calculated for the colli¬ 
sionless simulation because there was no star formation included 
and hence no new stars. 


Model 

Re 

Re, OS 

Re,ns 

cr 

0-5 


(kpc) 

(kpc) 

(kpc) 

(km s“^) 

(km s“^) 

nOlow 

7.4 

8.1 

6.1 

204 

212 

nOmed 

3.8 

5.7 

1.8 

263 

243 

nOhigh 

4.0 

6.1 

0.8 

226 

213 

nllow 

6.8 

6.8 

6.9 

224 

231 

nlmed 

3.8 

5.8 

1.8 

263 

242 

nlhigh 

4.2 

5.9 

1.9 

196 

185 

n2low 

5.9 

6.0 

5.8 

220 

224 

n2med 

4.0 

6.0 

2.0 

232 

205 

n2high 

5.0 

6.2 

2.3 

184 

176 

collisionless 

- 

7.2 

- 

155 

155 


I 2 OOII) . In a more pessimistic vein, ISpringe 3 ll2005h has re¬ 
cently shown that time integration methods commonly em¬ 
ployed in N-body simulation can lead to significant errors 
when using adaptive timesteps and highly eccentric orbits 
(see his Fig. 5). Thus, the central densities may be the re¬ 
sult of numerical errors which could be remedied by using 
new codes such as Springel’s GADGET-2. A much more 
thorough study of the time evolution of the surface density 
profile, the robustness to numerical integration parameters, 
and the dependence on various initial conditions will be re¬ 
quired to determine if merger remnants evolve into “regular” 
elliptical galaxies, and if so, which types of ellipticals. 

As a final comment, we remark that the centrally con¬ 
centrated new stellar population present in all models has 
affected even the collisionless component. The contraction 
of the collisionless component is clearly demonstrated in the 
lower-right plot of Figure 1101 The total surface density for 
n2high intersects the collisionless profile at R « 4.5 kpc, 
signaling that the gaseous dissipation and star formation 
included in n2high has contracted the mass profile in com¬ 
parison to the collisionless case. 


4■4-2 The Fundamental Plane 

To begin, we note that the previous section demonstrated 
that the collisionless stellar distribution responds to gaseous 
dissipation and star formation in the galaxy center by be¬ 
coming more concentrated. To quantify the contraction of 
stellar mass, we measure the half-mass radius Re of the total, 
the new, and the old stellar (bulge plus disk) mass and list 
these in Table [4] The half-mass radius is determined as the 
radius which encloses half the projected stellar mass. We do 
not attempt to compute a luminosity profile, nor do we try 
to fit a r^/'*-type luminosity profile since it does not appear 


to accurately characterize the surface density. Hence, Re is 
not identical to the effective radius, but it serves as a useful 
proxy. The half-mass radii of the new stars Re, ns shows that 
half of all stars formed during the medium and high feedback 
merger simulations end up within ~ 2 kpc of the center of 
the merger remnant. Further, Re, ns depends on n, the feed¬ 
back reservoir equation of state index, such that increasing 
n increases the half-mass radius of the newly formed stellar 
population. The dissipation, gas infall and star formation 
has affected the old stars as well. The half-mass radius of 
the old stars Re,os is ~ 14 — 20% smaller in the runs includ¬ 
ing gas, cooling, star formation, and feedback compared to 
the purely collisionless simulation. 

The central concentrations of newly formed stars affect 
the kinematics of the entire stellar population. To quantify 
the stellar kinematics, we calculate the one-dimensional ve¬ 
locity dispersion within an aperture centered on the poten¬ 
tial minimum. Table U lists the velocity dispersion within 
apertures of radius one-half Re and 5 kpc, averaged over a 
thousand projections. The typical scatter due to projection 
effects is ~ 10%. Because the one-dimensional velocity dis¬ 
persion is a measure of the depth of the potential well, the 
remnants with a smaller effective radius naturally lead to a 
larger velocity dispersion, and in all cases the velocity dis¬ 
persions are much larger than the collisionless simulation - 
the largest by over 50%. The dispersion profile also changes. 
In the collisionless case, the velocity dispersion is flat, i.e., 
the dispersion is identical when measured within an aperture 
of 3.6 or 5 kpc. In contrast, the medium and high feedback 
simulations have a velocity dispersion that increases toward 
the remnant center. 

While we reserve a more extensive study of these merger 
remnants to future work, it is interesting that the effective 
radii and velocity dispersions presented in Table |4] are con¬ 
sistent with the early-type Re-(J scaling relation found in 
the Sloan Digital Sky Survey (iBernardi et aLlF 20031 . see their 
Fig. 4 and 5). Further, the differences generated by the feed¬ 
back models scatter the merger remnants perpendicular to 
this scaling relation, suggesting that it is the interplay be¬ 
tween star formation and feedback which fixes the position 
of a merger remnant on the fundamental plane of hot stel¬ 
lar systems. The sensitivity of merger remnants to feedback 
models is potentially a useful discriminator for selecting the 
best model or parameters. 

Figure[TD]and Re,ns both demonstrate that the majority 
of new stars reside within the inner regions of the remnant. 
This is confirmed in Figure [TTJ which plots the cumulative 
mass distribution of new stars in each merger remnant. In all 
models approximately 80% of all new stellar mass is within 
10 kpc of the remnant center. In contrast, the simulations of 
MH96 led to merger remnants in which all of the new stellar 
material was within 1 kpc of the center. This suggests that 
our feedback model restricted the gas inflow and dispersed 
the central star-forming region. The new star mass cumu¬ 
lative profile is highly correlated with n. The n = 2 mod¬ 
els have the lowest cumulative mass distribution for both 
the high- and low-feedback series. Apparently, the stiffer 
equation of state has provided sufficient pressure support 
for stars to be formed over a larger range of radii, thus de¬ 
creasing the fraction of stars which reside within the inner 
~ 1 kpc of the remnant. 

Because stars are stochastically spawned from colli- 
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Figure 11. Cumulative new (formed during the merger simula¬ 
tion) stellar mass in the merger remnant for each model. Shown 
as solid, dotted and dashed thick lines that represent the n = 2, 
n = 1, and n = 0 feedback models, respectively. Then, the degree 
of feedback, i.e., high, med, and low, is denoted by the symbols 
as shown in the key. As in Figure 1101 the arrow indicates the 
gravitational softening length of 0.1 kpc. 


sional gas (see the last paragraphs in § I2.2II . it is possible 
that the new stellar component spuriously loses energy and 
angular momentum due to its artihcial coupling to the gas. 
This was explicitly checked by rerunning the n2med simu¬ 
lation with Ng = 10, i.e., each gas particle now spawns 10 
individual stellar particles as opposed to the fiducial value 
of 2. The resulting surface-density profile was very similar 
to n2med. The new stars were slightly more extended, with 
a half-mass radius Re,ns = 2.2 kpc, about 10% larger than 
n2med. The more extended new stars also affected the veloc¬ 
ity field and the dispersion was lower at half Re and higher 
at 5 kpc, but both of these effects were only ~ 3% of the 
values listed in Table 01 


5 DISCUSSION 

5.1 Convergence Study 

The resolution of N-body simulations such as those per¬ 
formed here is determined by the particle number N and the 
gravitational softening length h. Increasing N and decreas¬ 
ing h both result in higher resolution. Of course, this higher 
resolution comes at the expense of computational time. In 
practice, the selection of appropriate values for N and h is 
a balance between the resolution necessary to address the 
relevant processes and the available resources. In this sec¬ 
tion we explore variations in particle number and smoothing 
length for our fiducial merger simulation. SOO performed a 
similar study for the isolated disk. To isolate the effect of 
particle number, we select only one {n2med) of the feedback 
parameters to resimulate with the total particle number in¬ 
creased by factors of 2, 4 and 10. The ratio of dark-matter 
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Figure 12. Star-formation rate for our fiducial major merger 
with runs of varying resolution. 


to collisionless disk to collisional particle number remained 
fixed. The highest resolution simulation contained 1,000,000 
dark matter particles, 300,000 each of stellar disk and colli¬ 
sional gas, and 100,000 bulge particles for each galaxy, mak¬ 
ing the total number of particles used in this simulation 3.4 
million. The individual particle masses were 7.1 x 10® Mq 
for the dark matter, 1.8 x 10 ® Mq for the collisional gas, 
1.3 X 10® Mq for the stellar disk and 1.0 x 10® Mq for the 
stellar bulge. 

In general, increasing N allows the gravitational soften¬ 
ing length to be decreased as roughly However, be¬ 

cause of the increased computation time resulting from this 
scaling, we have chosen to use a smaller softening length for 
only one simulation. We chose our highest resolution simu¬ 
lation that had ten times the particle number. In this case, 
the softening length was 10^'^® (0.46) times smaller than the 
fiducial 0.1 kpc. 

Figure ITT] shows the star-formation rate for simulations 
of varied resolution. As the resolution is increased, either via 
N or h, we find a slight trend of lower star-formation rates 
during the initial passage, but the effects are less than 5% 
and, overall, there is excellent agreement between the global 
gas consumption. The consistency of star formation for all 
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resolution runs suggests that our fiducial particle number is 
sufficient to resolve the star-formation history. 

The merger remnants are also insensitive to the in¬ 
creased resolution. The surface density, half-mass radii, and 
velocity dispersions are all within 2% of the n2med simula¬ 
tion, with no trend as a function of higher resolution. How¬ 
ever, there are subtle differences between the high- and low- 
resolution simulations. In particular, the number of small 
condensations orbiting the central gas disk in the merger 
remnant (see bottom-right panel in Fig. 0 is increased by 
~ 30% in the higher resolution simulation. 


5.2 Comparison to Previous Work 

The simulations presented here extend prior studies of star 
formation and feedback in equal-mass gas-rich disk-galaxy 
major mergers. As in the work of MH96 and SOO, we find 
that mergers drive gas inflows which efficiently form stars 
and lead to a central concentration of newly formed stel¬ 
lar mass. Furthermore, we have quantified how these results 
depend on the assumed feedback efficiency and density de¬ 
pendence. However, the comparison to previous work is com¬ 
plicated by the many differences between the methodology 
used here and that employed in previous work. For clarity 
and definiteness, we find it useful to enumerate the differ¬ 
ences between the simulations presented here and those of 
MH96 and SOO: 

(i) To represent stellar feedback, MH96 imparted momen¬ 
tum “kicks” while SOO and this work chose to pressurize the 
ISM. In a recent paper. [Springel et al.l (l2005bh showed these 
“kicks” to be a weak form of feedback that was dependent 
upon the numerical resolution. We do not consider this fur¬ 
ther. 

(ii) Each simulation used a different disk-galaxy model. 

(iii) The star-formation recipe employed is identical in all 
simulations, but the normalization, i.e. the choice of c*, is 
different. The fiducial value used in all simulations here was 
0.03, while SOO and MH96 both used 0.004. For complete¬ 
ness, we will investigate a range of values both an order of 
magnitude smaller and larger than ours. 

(iv) MH96 represented the ISM as an isothermal gas, at 
10^ K. They argued that gas in disk galaxies efficiently cools 
and thus fixing the gas to this temperature scarcely affects 
the gas dynamics. SOO included adiabatic gas processes and 
shock heating as we do here. 

(v) Our simulations use a version of SPH which integrates 
the entropy, while MH96 and SOO used a version of SPH 
which integrated the energy. 

In order to better understand these differences, we will 
use this section to describe two types of tests. First, we will 
build the exact initial conditions of SOO and attempt to du¬ 
plicate his star-formation history. In addition, we will run 
his initial conditions with one of our models. In the second 
series, we use our fiducial model in a number tests which 
go through the above items, in order, and resimulate our 
fiducial n2med model with slight modifications aimed at il¬ 
lustrating the outcomes of each assumption. As mentioned 
above, we do not fully address item (i) because the compari¬ 
son to the “kicks” feedback used in MH96 has been rec ently 
performed by other authors (see lSpringel et ^l2005bll . We 



Time (/i'' Gyr) 

Figure 13. Star formation for model If of SOO. Shown are four 
runs, two that have the identical parameters as in SOO, which we 
have labeled accordingly, and the final two runs use the n2low 
and n2med parameters used throughout this work. Because the 
original work of SOO used the “energy, standard” formulation of 
SPH, our primary comparison run uses this as well. The other 
two runs use the newer version which is the standard for the 
main body of this paper. 

also do not perform a full exploration of the second differ¬ 
ence listed above, the disk galaxy model. A complete un¬ 
derstanding of the dependencies associated with the initial 
conditions is outside the scope o f this work (the beginning 
of such a project was broached in ICoxI ||2004|H . However, our 
first type of test, the explicit comparison using the models 
of SOO, will suffice to provide a basic understanding of item 
(ii) above. 

5.2.1 Comparison to SOO 

As a first task we’ll verify that our model can reproduce one 
of the quiescent star-formation histories and one of the merg¬ 
ing disk-galaxy star-formation histories previously shown in 
SOO. In order to perform this test, we must contruct identical 
initial disk models to those used in SOO. As a byproduct of 
this test, we will gain some insight as to how the disk galaxy 
model determines the resulting star-formation history. 

For the isolated galaxy comparison, we choose the II 
galaxy model outlined in § 4.1 of SOO. This model has V 200 = 
160 km s”*^, and thus a virial mass of 1.4 xlO^^ Mq. The 
halo has a concentration c = 5 and a spin parameter A = 
0.05. All of the baryons are in a disk, i.e., there is no bulge 
component, and 20% of the baryons are in a gaseous form. 
This results in mass fractions of md = 0.04, trig = 0.01, and 
fd = 0.2. Further, the gas and stellar disks are assumed to 
follow the same exponential distribution (i.e., a = 1) with 
a radial disk scale length of 6.4 kpc. The model is realized 
with 30,000 dark matter, 20,000 gas and 20,000 stellar disk 
particles for the same resolution used in SOO. 

When we introduced the parameter sets in § [23] we 
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provided the conversion to the exact model used in SOO. 
This parameter set had a lower star formation normaliza¬ 
tion c* = 0.004, but was otherwise very similar to our 
n2low. Although this wasn’t used for any of the models 
thus far, we now use it to simulate the galaxy model II 
in isolation. It is also significant to point out that the origi¬ 
nal SOO simulations used the “energy, standard” formulation 
of SPH where we hav e used a newer version introduced by 
ISoringel fc HernauistI (|2002h . We will simulate the SOO pa¬ 
rameter set with both here. As a reference to the parameter 
sets we’ve developed in this work, we will also run the II 
model with the nSlow and n2med parameter sets. 

Figure [13] shows the resulting star formation for the 
II galaxy model evolved 1.4 Gyr. As a confirmation of our 
ability to replicate the initial disk model, we note that both 
of the simulations that used the SOO parameter set have star- 
formation histories that are nearly identical to that shown in 
SOO (see his Fig. 7). Since the only difference between these 
two runs is the version of SPH, we conclude that the version 
of SPH makes little difference for evolving disk galaxies in 
isolation. 

The two simulations which evolved II with the n2low or 
n2med parameter set are quite different. The n2low model 
has a much higher star formation rate because it is nearly 
identical to SOO, only it has a larger value of c*. Because of 
this, nearly three times as much gas is converted to stars 
during the simulation and this rapid star formation quickly 
depletes the disk of high-density gas resulting in an sub- 
sequential decline of the star-formation rate. The average 
star-formation rate of the simulation that uses the n2med 
parameter set is much closer to that of the SOO runs, how¬ 
ever there are strong oscillations present. These oscillations 
were seen in the med- and high-feedback runs during the 
Sbc major merger and were discussed briefly at the end of 

§0 

An interesting comparison can be made between the 
star formation reported in Figure[T3]for the isolated II disk, 
and that shown earlier in Figure[T]for the isolated Sbc disk. 
First, with the n2lo'w parameter set, we note that the Sbc 
disk has a steady star-formation rate of ~ 0.75 M 0 yr“^. 
This same parameter set for the II disk begins with a much 
larger star-formation rate (~ 6 M 0 yr“^) that quickly de¬ 
clines as the high-density gas is consumed. Not surprisingly, 
one concludes that the star-formation history is a direct con¬ 
sequence of the amount and distribution of high-density gas, 
bot h of which are a byproduct of the initial d isk model (see 
also lLi et al.ll200Rl2006l : [Springel et al.ll2005b[l . This leads us 
to suggest a natural extension to our work: extend our tests 
to a large set of initial disk models with properties drawn 
from observations, as with our Sbc-disk model. By compar¬ 
ing these tests to observed relations between Hubble-type 
and star formation, one may provide further constraints on 
the star formation and feedback model. 

After assuring ourselves that we are able to evolve the 
exact II model of SOO in a similar fashion and investigating 
the differences induced by either the version of SPH, or the 
alternate star formation and feedback model, we now move 
on to test one of his major mergers. For this test, we select 
the bulgeless A1 major merger. These initial disk models are 
slightly smaller than the II disk, with V 200 = 120 km s“^, 
and thus their mass is 6.0 xlO^^ M 0 . The halo has a con¬ 
centration c = 5 and a spin parameter A = 0.05. As with 



Time (/i'' Gyr) 

Figure 14. Star formation for the A1 bulgeless major merger 
of SOO. As with Figure 1131 we show one run that is designed 
to reproduce the star-formation history of SOO, and three other 
runs which have small variants in either the star formation and 
feedback parameters, or the version of SPH. 


the II model, all of the baryons are in a disk, i.e., there is no 
bulge component, and the mass fractions are = 0.04 and 
rrig = 0 . 01 , where the gaseous and stellar disks are similarly 
distributed. This model is realized with an identical particle 
number as II. 

Figure [TT] shows the star formation rate for the A 1 ma¬ 
jor merger. As with the II test, we run four simulations: 
one identical to SOO, one with the same star formation and 
feedback parameters as SOO but with the newer version of 
SPH, and two with our parameter sets n2low and n2med. 
Once again, when we use the SOO parameter set and the 
identical version of SPH, we can reproduce the A1 merger 
star-formation rate exactly (see Fig. 14 of SOO). 

In contrast to the isolated case, when we use the 
SOO parameter set with the newer version of SPH, the 
star-formation history shows marked changes. The star- 
formation rate in this case is oscillatory and suppressed rela¬ 
tive to the SOO result. These changes occur primarily during 
the bursts of star formation that result after the hrst passage 
and the final merger and are likely a result of the improved 
entropy conservation in the newer version of SPH. We will 
perform more tests along this line, and discuss this further 
in §[523] 

The other two runs, which use the n2low and n2med 
parameter sets, also show distinct differences from SOO. As 
in the case of the isolated disk, n2low shows much larger star 
formation during the quiescent stages of evolution. However, 
the bursts are relatively similar in amplitude, especially to 
the SOO model with the same SPH version. The n2med run 
shows the oscillations that were also present in II and Sbc, 
and there is almost no burst of star formation subsequent to 
the first passage at T = 1.0— 1.7 h~^ Gyr. The hnal merger, 
however, contains the highest star-formation rate of any of 
the four models. 
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Figure 15. Star-formation rate during the fiducial merger for 
different values of c*, the parameter that sets the star formation 
efficiency. 

From the tests performed in this section, we arrive at 
three primary conclusions. First and foremost, our model 
reduces to that of SOO, and as such, if we use his initial 
conditions, we can reproduce his results. Second, we note 
that the increased value for the star-formation normaliza¬ 
tion c*, unsursprisingly, produces an increased star-formtion 
rate. This increase occurs only during the quiescent evolu¬ 
tion. During the merger-induced bursts, the star-formation 
rate is nearly equivalent, as demonstrated in Figure [TT] The 
third conclusion pertains to the version of SPH. While this 
made very little difference while evolving the quiescent disk, 
there are substantial differences that occur during the more 
dynamic merger. We will perform a more thorough explo¬ 
ration of these last two dependencies, c* and the version of 
SPH, using our Sbc disk-galaxy model in the sections that 
follow. 


5.2.2 Star-Formation Recipe: Normalization 

The parameter c* normalizes the course-grained recipe to 
convert gas to stars via Equation ([T}. While H4.ll showed 
that c* = 0.03 reproduced the observational Kennicutt law 
m over two orders of magnitude in gas surface density, the 


Table 5. Star-formation properties of major merger simulations 
as a function of c*. Quantities are defined in Table [3] 


C* 

e 

SFRmax 

(Mgyr-i) 

Tlirg 

(Gyr) 

0.3 

0.51 

82 

1.30 

0.06 

0.49 

57 

1.36 

0.03 (n2med ) 

0.47 

56 

1.34 

0.015 

0.45 

38 

1.21 

0.004 

0.28 

22 

0.01 


observational scatter appears large enough to allow for a 
range of star-formation efficiencies. In this section, we res¬ 
imulate the fiducial n2med model with c*= 0.004, 0.015, 
0.06, and 0.3. These values of c* represent excursions from 
our fiducial value (0.03) of an order of magnitude (0.004 and 
0.3) and a factor of two (0.015 and 0.06). Both MH96 and 
SOO used c* « 0.004. 

Figure [15] shows the star-formation rate for the fidu¬ 
cial n2med simulation and the four mergers with alternate 
star-formation normalizations. Not surprisingly, the star for¬ 
mation prior to the final merger correlates with the star 
formation efficiency c*, but the situation reverses after the 
merger. In fact, the low efficiency of star formation means 
very little gas is consumed prior to the final merger, and by 
the end of the simulation, at T = 3.0 Gyr, the c* = 0.004 
simulation has the highest rate of star formation. In a man¬ 
ner similar to H4.1I the star formation can be quantified by 
the gas consumption, the maximum star-formation rate, and 
the duration in which the merger would qualify as a LfRG, 
which are listed in Table [5| 

Varying the star-formation efficiency highlights the dif¬ 
fering effects of star formation versus feedback. For exam¬ 
ple, changing the fiducial n2med model by decreasing c* to 
0.004 (as done above) or increasing rfb to 2.0 {n2high), re¬ 
duces SFRmax and e. Increasing the feedback, as in n2high, 
restricts the ability of gas to reach high densities, but gas 
which does so will be rapidly converted to stars. Hence, the 
gas consumption e is strongly affected while the peak star- 
formation rate is only slightly less than n2med. In contrast, 
decreasing c* only slightly alters the gas density but strongly 
reduces the star-formation efficiency of the high-density gas. 
In this case, the maximum star-formation rate is less than 
half n2med, but the gas consumption e is higher than when 
increasing the feedback. 

Because the star-formation rate was affected by differ¬ 
ent values of the star-formation efficiency c*, it’s reasonable 
to ask if the agreement with the Kennicutt law has been de¬ 
stroyed. Figure [16] answers this question by showing the gas 
and star-formation rate surface densities in a 2 kpc aperture 
for each of the simulations with different c*. As expected, 
due to our coupling of star formation to the gas density via 
Equation 0, the slope of each simulation tracks the ob¬ 
served Kennicutt law. However, the various normalizations 
have shifted the amplitude of star formation so that the val¬ 
ues 0.3, and 0.004 are on the edges of the observed range. 


© 0000 RAS, MNRAS 000, 000-000 





















Feedback in simulations of disc-galaxy major mergers 25 




Log S gas (Mo pc"') Log 

^ gas (Mo pc"') 


Figure 16. Identical to Figures [8] except shown here are 
simulations with varying star-formation normalization c*. The 
(red) squares are the fiducial n2med simulation and were sim¬ 
ilarly shown in Figure O The (blue) circle, (orange) aster¬ 
isk, (yellow) cross, and (black) triangle symbols are for c* = 
0.3,0.06,0.015, and 0.004, respectively. 


Figure 17. This figure is identical to Figure [3] except here we 
compare our fiducial n2low model to variants that have an altered 
star-formation law according to Equation II23II . with N = Q, 2 
(opposed to the fiducial N = 1.5). All simulations are run with 
the threshold density for star-formation set to zero. 


5.2.3 Star-Formation Recipe: Slope 


The previous section demonstrated that the normalization 
of the pro jected star - forma tion law, and hence the agree¬ 
ment with KennicuttI (Il998l ) , is set by the free parameter c*. 
However, we have still not demonstrated what sets the slope 
of the star-formation law. Previous efforts have argued that 
the observed star-formation law is a natural conseq uence of 
the turbulent interstellar medium llKravtsovll2003l ) or that 
it res ults from gravi tational instabilities in gaseous spiral 
disks llLi et al.ll2004 l. In order to investigate the origin of 
the slope of the star-formation law, we will follow a pro¬ 
cedure similar to iMihos et al.l (Il99ll) and modify our star- 
formation prescription, i.e.. Equation ©• To this end, we 
write a generalized star-formation law as 


dp* ^ ^ N 
At 


(23) 


where C is a dimensionful free parameter and N sets the 
gas-density power-law dependence. In our fiducial model, 
C = VdyrGc*, and N = 1.5. As variants, we will try N = Q 
and 2, and fix C such that the global star-formation rate is 
identical at pth- To simplify the comparison, we only show 
results from simulations where the star formation is not fixed 
to be zero below pth- We also performed identical simulations 
that included the threshold density, however in this case the 
star-formation surface density quickly plunged to zero below 
~ 10 Mq pc“^, as in Figure |3] and thus there was very little 
difference between the various star-formation laws. 

Figure fTTl shows the star-formation law for our isolated 
Sbc disk galaxy for the three variants of the star formation 
prescription. We see a clear trend that the slope of the pro¬ 
jected star-formation law follows the power of N in Equa- 
tion (1231) . This correlation was also noted by iMihos et al.l 


(119911) . who additionally showed that the steeper powers of 
N result in starbursts of greater intensity. From these con¬ 
siderations we conclude that the star-formation law, both in 
normalization and slope, is fixed by the empircal Kennicutt 
law. 


5 . 2.4 Isothermal ISM 

Next we consider differences generated by assuming the ISM 
is an isothermal gas with temperature 10“* K. Because of the 
efficient cooling in quiescent disk galaxies, it is reasonable 
to assume their thermal energy is maintained near 10^ K. 
Indeed, this is precisely the case for our model galaxy when 
simulated in isolation, there is no difference. 

However, galaxy mergers are a much more violent event 
in which there i s ample shock heating of low-density gas 
(ICox et al.l I 2 OO 4 I ) . This was demonstrated in m\ where 
phase diagrams showed large amounts of hot (~ 10® K), 
low density (10~^ — 1O“^M0 pc“®) gas. In this state, the 
cooling time is much longer than the orbital time, and the 
gas becomes locked in this phase, unable to cool. If the tem¬ 
perature of this gas were fixed at 10"^ K, this gas would be 
much more likely to fall back to the galactic center and could 
significantly affect the star-formation rate. 

To investigate this possibility, we have run our n2med 
model with the thermal gas temperature (the u reservoir) 
fixed to 10"^ K. (In order to maintain a stable isolated disk, 
the feedback reservoir q still operates normally.) Under these 
assumptions, the star-forming gas is scarcely different than 
n2med while all non-star forming gas is effectively 10“^ K, 
without gas which has recently formed stars and hasn’t ther- 
malized its turbulent energy yet. The resulting star forma¬ 
tion, in comparison to n2med, is shown in Figure 1181 The 
isothermal assumption results in a maximum star-formation 
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Figure 18. Star-formation rate for our n2med model and for the 
identical model except the temperature of non-star forming gas 
is fixed at 10^ K. 

rate that is four times that of n2med and increases the 
gas consumption e by 20% over n2med. This gas fraction 
is nearly identical to the quantity of gas which ends up in 
the hot, low-density phase in n2med. Further, the periods of 
maximum difference between the two star-formation rates is 
precisely when large shock heating occurs. It thus appears 
that assuming the ISM is isothermal at 10^ K grossly under¬ 
predicts the amount of shock heating due to the merger and 
results in much more cold gas available for star formation. 


5.2.5 SPH Version 

All simulations presented up to this point have used a ver¬ 
sion of SPH which integrates the entropy and was formu- 
late d to conserve both energy a nd entropy, when appropri¬ 
ate dSpringel fc Hernauisdl^02li . In order to discern the ef¬ 
fects of our choice of SPH we resimulate the fiducial n2med 
merger with alternate versions. Specifically, we use the en¬ 
ergy integrating “energy, geometric”, “energy, standard”, 
and “energy, asymmetric” versions of SPH, as well as the en¬ 
tropy integrating “entropy, sta ndard” version, where we fol - 
low the naming convention of ISpringel fc Hernauisd (l2002l l 
and refer the interested reader to their paper for further de- 
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Figure 19. A comparison of our fiducial n2med model run with 
different versions of SPH commonly used in numerical simula¬ 
tions. 


Table 6. Star formation properties of major merger simulations 
as a function of SPH formulation. Quantities are defined in Ta¬ 
ble |3] 


SPH Version 

e 

SFRmax 

(M©yr-i) 

Tlirg 

(Gyr) 

entropy (n2med) 

0.47 

56 

1.34 

energy, standard 

0.58 

55 

1.48 

energy, geometric 

0.65 

65 

1.52 

energy, asymmetric 

0.58 

61 

1.42 

entropy, standard 

0.52 

60 

1.38 


tails. The galaxy-merger simulations of SOO and MH96 used 
the “energy, standard” version. 

Figure [19] shows the star-formation rate during the 
SPH runs. There are two notable features present in the 
star-formation histories owing to the version of SPH. The 
first feature is the systematic enhancement of star forma¬ 
tion from T = 0.6 to 1.6 Gyr i n runs that do not use 
the fiducial ISpringel fc HernauistI (|2002| ~) version of SPH. 
The second feature is present after the final merger, from 
T = 1.8 to 2.0 Gyr, when the energy versions of SPH lead 
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Table 7. Star formation properties of major merger simulations 
as a function of ISM metallicity. The metallicity is given relative 
to solar and the remaining quantities are defined in Table l3l 


Metallicity 

(Z©) 

e 

SFRmax 

(Moyr-i) 

Tlirg 

(Gyr) 

zero {n2med) 

0.47 

56 

1.34 

10-3 

0.49 

54 

1.38 

10-3 

0.50 

54 

1.40 

10-1.5 

0.54 

47 

1.44 

10-1 

0.58 

67 

1.44 

1 

0.60 

53 

1.52 


to a large peak of star formation while the entropy versions 
{n2med and “entropy, standard”) of SPH have an oscilla¬ 
tory star-formation rate. In general, the entropy versions 
lead to lower star-formation rates and less gas consumption. 
We speculate that these features arise because the entropy 
versions of SPH more accurately treat point-like energy in¬ 
jections and adiabatic heating as gas falls into the central 
region of the halo, as demonstrated in ISpringel fc HernauistI 
(1200211 ■ and therefore may be more accurate than energy ver¬ 
sions of SPH. Additionally, while versions of SPH not tested 
here may produce still different results, we feel that the tests 
performed here survey the plausible range of outcomes and 
accurately reflect the modeling uncertainties resulting from 
versions SPH. 

While these tests show that the version of SPH can 
alter the star formation during a galaxy interaction, we note 
that these differences are smaller than differences induced 
by other modeling assumptions. For instance, Table|6]shows 
that the total amount of gas consumed and maximum star- 
formation rate for this series of mergers varies by < 20%. 
By comparison, changes in the star formation efficiency (c*) 
or feedback parameters (rfb and n) can produce differences 
of over a factor of 2 (see Tables 13 and O . 


5.3 ISM Metallicity 

The ISM radiatively cools at a rate which is highly depen¬ 
dent upon the metallicity. Up to now, cooling has been cal¬ 
culated assuming the ISM is a primordial plasma with zero 
metallicity even though metals are tracked by the code, as 
mentioned in §[I1 In this section we resimulate the fidu¬ 
cial n2med model with rad iative cooling calculated from the 
tabulated cooling curves of ISutherland fc Dopital (Il993l l and 
assumed gas metallicities, in units of solar, of 10“®, 10“^, 
10“^'®, and 1. As with our fiducial model, we continue 
to ignore spatial and temporal changes in the metallicity. 
The resulting star-formation rate is shown in Figure [20] As 
expected, the increased metallicity leads to a higher fraction 
of the gas in the cold phase, and increases the star-formation 
rate. The more efficient metal-line cooling is strongest for gas 
which is ~ 10® K, and hence the star formation differences 
are not apparent until gas begins to populate this temper¬ 
ature range. Gas at these temperatures is produced by the 
interaction itself, in shocks that occur between gas in the 
progenitor disks. While the maximum star-formation rate 
does not dramatically change for simulations with increased 
gas metallicity, the gas consumption e steadily increases. In 
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Figure 20. Star-formation rate for our fiducial major merger 
with various assumptions for the inter stellar medium metallicity. 
The metallicities listed are in solar units. 


fact, there is a 13% increase in e during the merger with 
solar versus zero metallicity. The star-formation properties 
are listed in Table [7| 

The variations in gas consumption suggest that appre¬ 
ciable errors can be generated when numerical simulations 
treat the ISM as a primordial plasma. While these varia¬ 
tions are smaller than those induced by altering the star 
formation and feedback parameters (see Tables [3] and [S]l 
the cooling rate can be calculated fully self-consistently be¬ 
cause the simulation tracks metal enrichment (see 3231). 
Toward this end a number of groups have recently in¬ 
cluded metallicity-dependent cooling based upon a model 
for t he enrichment of metals from supernovae (SN ) (see , 
e.g.. ISca nnapiec o et al.l l2 005l. l2006l: lOka moto et al] l2005l : 


iTornatore et al.l 2004 ; Mosconi et al l 12001 ) and have found 
similar trends to those reported here. Unfortunately, these 
models necessarily require a number of assumptions regard¬ 
ing the time delay of SNIa, the physical scale that stellar 
particles distribute their metals, the diffusion of metals, and 
the scale in which metals are smoothed when the cooling is 
calculated. At the current moment it is unclear how best to 
implement the physics of metal enrichment and metallicity- 
dependent cooling but these physical processes are certain to 
garner a significant amount of attention in the near future. 
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6 CONCLUSION 

We have carried out a series of isolated and merging disk 
galaxy simulations which include simple recipes for imple¬ 
menting star formation and feedback. The free parameters 
of these models are selected based upon the stability of our 
isolated disk galaxy model and excursions are explored with 
increased and descreased stability. These models were then 
used during equal-mass prograde major mergers to inves¬ 
tigate the resulting star-formation and remnant properties. 
The primary conclusions can be summarized as follows: 

• The slope and normalization of the star-formation law 
are fixed by adopting a recipe in which stars are formed 
from high-density gas according to p* « c*Pgas. The slope 
is a result of the 1.5 index, and the normalization is set by 
c*. Star formation designed in this manner, independent of 
the feedback model, can naturally reproduce the empirical 
Kennicutt law. The above holds for both quiescent and in¬ 
teracting galaxies. 

• Provided that star formation is designed as in the first 
item above, the rate at which stars are formed is then de¬ 
termined by the feedback efficiency rfb. The star-formation 
rate in isolated disk galaxy simulations is largely insensi¬ 
tive to the equation of state at high densities (in our model 
this is controlled by n), because any feedback sufficient to 
stabilize the galactic disk is also efficient enough to confine 
star-forming gas to a narrow range of densities around the 
threshold density pth. 

• During a merger, where gas can acheive much higher 
densities than in the isolated case, the star formation rate 
depends on both rfb and n. However, the overall gas con¬ 
sumption is fixed primarily by rn,, the feedback efficiency. 
Because rfb determines the pressurization of star-forming re¬ 
gions around the threshold density, large values can prevent 
gas from ever reaching star-forming densities. 

• The maximum star-formation rate during the merger 
simulations depends on rn, in a way which is reasonably 
well approximated by a toy model of star-forming gas in a 
fixed potential (Equation I22II . 

Our results extend, and are consistent with, previous 
numerical simulations of star formation in colliding disk 
galaxies. A detailed comparison finds that previous versions 
of SPH predict an increased gas consumption by varying 
amounts. We also determined that an isothermal gas as¬ 
sumption at 10“* K also predicts much more efficient star 
formation during a major merger. 

The above items collectively imply that requiring free 
parameters to yield a stable isolated disk and star forma¬ 
tion that is consistent with the Kennicutt law provides 
a fairly tight constraint on the star-formation recipe but 
does not uniquely determine the feedback parameters. The 
large amount of freedom in selecting the feedback param¬ 
eters makes a significant difference to the maximum star- 
formation rate during a galaxy merger, while still fulfilling 
these basic requirements. 

Generally speaking, the integrated gas consumption 
during a merger is less sensitive to the details of the feed¬ 
back parameters, due to a trade-off between star-formation 
rate and the duration of the burst. In models with highly 
efficient feedback, strong bursts of star-formation are sup¬ 
pressed. However, star formation then continues until all the 


gas at star-forming densities has been consumed, leading to 
a similar overall gas consumption. 

While we argue that the formulation of star forma¬ 
tion in the numerical simulations is well defined, the de¬ 
termination of the feedback parameters is not and thus 
we ask the question: Is there any hope of determing a 
unique set of feedback parameters which may then be 
used for a large series of galaxy merger simulations? In 
our (optimistic) opinion, the answer to this question is 
yes. However, doing so may be a long, arduous, and ill- 
defined process. As a first step we note that there is a 
strong observationa l link between galaxy interactions and 
star formation(e.g.. iBor ne et al.ll200(]l: iLambas et al.l l2003l : 


iBarton Gillespie et ^ 20031 : Nikolic et al.l 2004 ). Theoret¬ 

ically, models of galaxy formation with merger induced 
star -formation naturally p roduce the Lyman-break galax- 


star -tormation naturally p roduce tne Lyman-break galax¬ 
ies dSomerville et al. 200l|) and possibly also submillimeter 
galaxies ( Baugh et al.l 120051') . Thus, it seems reasonable to 
require that our merging gas-rich spiral galaxies induce a 
significant burst of star formation. In this sense, the “high” 
feedback models are disfavored. Further constraints on burst 
strengths and ages may place tighter bounds on feedback 
models. For instance, the prolonged bursts present in the 
n = 2 models may be at odds with observations that sug¬ 
gest starburst have a shorter duration ~ 100 Myr. We also 
note that it would be very helpful to have a statistical obser¬ 
vational catalog of mergers and merger remnants, including 
estimates of the orbits and galaxy masses, burst strengths 
and ages, and any associated kinematics to compare with 
the simulations. 

Ultimately, one of the most promising tests of star for¬ 
mation and feedback models may require detailed compar¬ 
isons between observations of quiescent or interacting galax¬ 
ies and simulations designed to reproduce these galactic sys¬ 
tems. Modeling of specific galaxy systems has already bee n 
done in several instances such a s “The M ice” (lBarnes| [2004 ), 
M51 dSalo fc Laurikainenll200oh. ArpllS (lLamb_et_^lJl99S 


Arpl l9 ( Hearn fc Lamb 2001 ), Arp22 0 (iMcDowell et al 
2003), NGG 7714/ 15 ds truck fc SmithI [ 20031 : ISmith et al 


2005li NGC 2442 jMihos fc BothunI Il997ll. IC2163/NGG 


2207 dstruck et al.ll2005l l and NGG 7252 lllibbard fc MihosI 
Il995ll . although none of these explored a range of feedback 
parameterizations. It is evident that much work must be 
done to determine what initial conditions are plausible and 
which physical processes must be included to reproduce the 
varied and extreme environments of galaxy mergers. 

Even with the ambiguities present in our models it is 
encouraging that under a wide range of assumptions we are 
able to produce stable, gas-rich quiescent galaxies and large 
prolonged bursts of star formation. Furthermore, for the 
models that produce large bursts of star formation, the in¬ 
sensitivity of the global gas consumption to the details of 
the model leads us to believe that we can apply these tech¬ 
niques to a large series of minor as well as major merger 
simulations and achieve a reliable estimate of the ability of 
galaxy interactions to drive star formation, a project which 
has already been started. While the same universality does 
not hold for the maximum star-formation rate of each model, 
this is something that might allow us to constrain feedback 
by matching to observations, as we mentioned above. 

In addition, we have presented evidence that the rem¬ 
nant properties, such as half-mass radius and velocity disper- 
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sion, are highly dependent upon the star formation and feed¬ 
back assumptions. Low feedback models are more concen¬ 
trated and have higher central velocity dispersions. Hence, 
it is likely that performing a large series of runs, and com¬ 
paring them to a statistical sample of elliptical galaxies, will 
significantly discriminate between the models. It is also pos¬ 
sible that some leverage can be garnered from the signihcant 
differences in gas properties predicted by the different feed¬ 
back models. By comparing phase diagrams such as those in 
i)4.3l to mass-weighted phase diagrams of observed galaxies 
there may be some models that could be ruled out. Another 
avenue for comparison between the models and observations 
is the gas morphology. For instance, the feedback parame¬ 
ters affect the vertical structure of the gas disk, and hence 
quan tities like the inclination dependence of dust attenua¬ 
tion lljonsso nl l2004ll . 
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